perfectGasI.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-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 "perfectGas.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Specie>
31 inline Foam::perfectGas<Specie>::perfectGas(const Specie& sp)
32 :
33  Specie(sp)
34 {}
35 
36 
37 template<class Specie>
39 (
40  const word& name,
41  const perfectGas<Specie>& pg
42 )
43 :
44  Specie(name, pg)
45 {}
46 
47 
48 template<class Specie>
51 {
53 }
54 
55 
56 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
57 
58 template<class Specie>
59 inline Foam::scalar Foam::perfectGas<Specie>::rho(scalar p, scalar T) const
60 {
61  return p/(this->R()*T);
62 }
63 
64 
65 template<class Specie>
66 inline Foam::scalar Foam::perfectGas<Specie>::h(scalar p, scalar T) const
67 {
68  return 0;
69 }
70 
71 
72 template<class Specie>
73 inline Foam::scalar Foam::perfectGas<Specie>::Cp(scalar p, scalar T) const
74 {
75  return 0;
76 }
77 
78 
79 template<class Specie>
80 inline Foam::scalar Foam::perfectGas<Specie>::e(scalar p, scalar T) const
81 {
82  return 0;
83 }
84 
85 
86 template<class Specie>
87 inline Foam::scalar Foam::perfectGas<Specie>::Cv(scalar p, scalar T) const
88 {
89  return 0;
90 }
91 
92 
93 template<class Specie>
94 inline Foam::scalar Foam::perfectGas<Specie>::sp(scalar p, scalar T) const
95 {
96  return -this->R()*log(p/Pstd);
97 }
98 
99 
100 template<class Specie>
101 inline Foam::scalar Foam::perfectGas<Specie>::sv(scalar p, scalar T) const
102 {
104  return 0;
105 }
106 
107 
108 template<class Specie>
109 inline Foam::scalar Foam::perfectGas<Specie>::psi(scalar p, scalar T) const
110 {
111  return 1.0/(this->R()*T);
112 }
113 
114 
115 template<class Specie>
116 inline Foam::scalar Foam::perfectGas<Specie>::Z(scalar p, scalar T) const
117 {
118  return 1;
119 }
120 
121 
122 template<class Specie>
123 inline Foam::scalar Foam::perfectGas<Specie>::CpMCv(scalar p, scalar T) const
124 {
125  return this->R();
126 }
127 
128 
129 template<class Specie>
130 inline Foam::scalar Foam::perfectGas<Specie>::alphav(scalar p, scalar T) const
131 {
132  return 1/T;
133 }
134 
135 
136 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
137 
138 template<class Specie>
140 {
141  Specie::operator+=(pg);
142 }
143 
144 
145 template<class Specie>
146 inline void Foam::perfectGas<Specie>::operator*=(const scalar s)
147 {
148  Specie::operator*=(s);
149 }
150 
151 
152 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
153 
154 template<class Specie>
155 inline Foam::perfectGas<Specie> Foam::operator+
156 (
157  const perfectGas<Specie>& pg1,
158  const perfectGas<Specie>& pg2
159 )
160 {
161  return perfectGas<Specie>
162  (
163  static_cast<const Specie&>(pg1) + static_cast<const Specie&>(pg2)
164  );
165 }
166 
167 
168 template<class Specie>
169 inline Foam::perfectGas<Specie> Foam::operator*
170 (
171  const scalar s,
172  const perfectGas<Specie>& pg
173 )
174 {
175  return perfectGas<Specie>(s*static_cast<const Specie&>(pg));
176 }
177 
178 
179 template<class Specie>
180 inline Foam::perfectGas<Specie> Foam::operator==
181 (
182  const perfectGas<Specie>& pg1,
183  const perfectGas<Specie>& pg2
184 )
185 {
186  return perfectGas<Specie>
187  (
188  static_cast<const Specie&>(pg1) == static_cast<const Specie&>(pg2)
189  );
190 }
191 
192 
193 // ************************************************************************* //
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Perfect gas equation of state:
Definition: perfectGas.H:96
scalar Cv(scalar p, scalar T) const
Return Cv contribution [J/(kg K].
Definition: perfectGasI.H:87
scalar psi(scalar p, scalar T) const
Return compressibility [s^2/m^2].
Definition: perfectGasI.H:109
scalar alphav(const scalar p, const scalar T) const
Return volumetric coefficient of thermal expansion [1/T].
Definition: perfectGasI.H:130
scalar e(const scalar p, const scalar T) const
Return internal energy contribution [J/kg].
Definition: perfectGasI.H:80
autoPtr< perfectGas > clone() const
Construct and return a clone.
Definition: perfectGasI.H:50
perfectGas(const Specie &sp)
Construct from components.
Definition: perfectGasI.H:31
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
Definition: perfectGasI.H:59
void operator+=(const perfectGas &)
Definition: perfectGasI.H:139
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
Definition: perfectGasI.H:123
scalar h(const scalar p, const scalar T) const
Return enthalpy contribution [J/kg].
Definition: perfectGasI.H:66
scalar Cp(scalar p, scalar T) const
Return Cp contribution [J/(kg K].
Definition: perfectGasI.H:73
scalar sv(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cv/T [J/kg/K].
Definition: perfectGasI.H:101
scalar sp(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cp/T [J/kg/K].
Definition: perfectGasI.H:94
scalar Z(scalar p, scalar T) const
Return compression factor [].
Definition: perfectGasI.H:116
void operator*=(const scalar)
Definition: perfectGasI.H:146
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.
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)
volScalarField & p