Loading...
Searching...
No Matches
perfectGasI.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
28#include "perfectGas.H"
30
31using namespace Foam::constant::thermodynamic;
32
33// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
34
35template<class Specie>
36inline Foam::perfectGas<Specie>::perfectGas(const Specie& sp)
37:
38 Specie(sp)
39{}
40
41
42// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43
44template<class Specie>
46(
47 const word& name,
48 const perfectGas<Specie>& pg
49)
51 Specie(name, pg)
52{}
53
54
55template<class Specie>
58{
60}
61
62
63template<class Specie>
65(
66 const dictionary& dict
67)
68{
73// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
74
75template<class Specie>
76inline Foam::scalar Foam::perfectGas<Specie>::rho(scalar p, scalar T) const
77{
78 return p/(this->R()*T);
79}
80
81
82template<class Specie>
83inline Foam::scalar Foam::perfectGas<Specie>::H(scalar p, scalar T) const
84{
85 return 0;
86}
87
88
89template<class Specie>
90inline Foam::scalar Foam::perfectGas<Specie>::Cp(scalar p, scalar T) const
91{
92 return 0;
93}
94
95
96template<class Specie>
97inline Foam::scalar Foam::perfectGas<Specie>::E(scalar p, scalar T) const
98{
99 return 0;
100}
101
102
103template<class Specie>
104inline Foam::scalar Foam::perfectGas<Specie>::Cv(scalar p, scalar T) const
105{
106 return 0;
107}
108
109
110template<class Specie>
111inline Foam::scalar Foam::perfectGas<Specie>::S(scalar p, scalar T) const
112{
113 return -this->R()*log(p/Pstd);
114}
115
116
117template<class Specie>
118inline Foam::scalar Foam::perfectGas<Specie>::psi(scalar p, scalar T) const
119{
120 return 1.0/(this->R()*T);
121}
122
123
124template<class Specie>
125inline Foam::scalar Foam::perfectGas<Specie>::Z(scalar p, scalar T) const
126{
127 return 1;
128}
129
130
131template<class Specie>
132inline Foam::scalar Foam::perfectGas<Specie>::CpMCv(scalar p, scalar T) const
133{
134 return this->R();
135}
136
137
138// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
139
140template<class Specie>
142{
143 Specie::operator+=(pg);
144}
145
146
147template<class Specie>
148inline void Foam::perfectGas<Specie>::operator*=(const scalar s)
149{
150 Specie::operator*=(s);
151}
152
153
154// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
155
156template<class Specie>
157inline Foam::perfectGas<Specie> Foam::operator+
158(
159 const perfectGas<Specie>& pg1,
160 const perfectGas<Specie>& pg2
161)
162{
163 return perfectGas<Specie>
164 (
165 static_cast<const Specie&>(pg1) + static_cast<const Specie&>(pg2)
166 );
167}
168
169
170template<class Specie>
171inline Foam::perfectGas<Specie> Foam::operator*
172(
173 const scalar s,
174 const perfectGas<Specie>& pg
175)
176{
177 return perfectGas<Specie>(s*static_cast<const Specie&>(pg));
178}
179
180
181template<class Specie>
182inline Foam::perfectGas<Specie> Foam::operator==
183(
184 const perfectGas<Specie>& pg1,
185 const perfectGas<Specie>& pg2
186)
187{
188 return perfectGas<Specie>
189 (
190 static_cast<const Specie&>(pg1) == static_cast<const Specie&>(pg2)
191 );
192}
193
194
195// ************************************************************************* //
#define R(A, B, C, D, E, F, K, M)
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition autoPtr.H:65
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition dictionary.H:133
Perfect gas equation of state.
Definition perfectGas.H:87
scalar Cv(scalar p, scalar T) const
Return Cv departure [J/(kg K].
Definition perfectGasI.H:97
static autoPtr< perfectGas > New(const dictionary &dict)
Definition perfectGasI.H:58
scalar E(const scalar p, const scalar T) const
Return internal energy departure [J/kg].
Definition perfectGasI.H:90
scalar psi(scalar p, scalar T) const
Return compressibility rho/p [s^2/m^2].
scalar H(const scalar p, const scalar T) const
Return enthalpy departure [J/kg].
Definition perfectGasI.H:76
autoPtr< perfectGas > clone() const
Construct and return a clone.
Definition perfectGasI.H:50
perfectGas(const Specie &sp)
Construct from components.
Definition perfectGasI.H:29
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
Definition perfectGasI.H:69
scalar S(const scalar p, const scalar T) const
Return entropy [J/(kg K)].
void operator+=(const perfectGas &)
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
scalar Cp(scalar p, scalar T) const
Return Cp departure [J/(kg K].
Definition perfectGasI.H:83
scalar Z(scalar p, scalar T) const
Return compression factor [].
void operator*=(const scalar)
A class for handling words, derived from Foam::string.
Definition word.H:66
volScalarField & p
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Thermodynamic scalar constants.
const scalar Pstd
Standard pressure: default in [Pa].
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
dimensionedScalar log(const dimensionedScalar &ds)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
Definition exprTraits.C:127
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dictionary dict