rhoTabulatedI.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) 2020-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 "rhoTabulated.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Specie>
32 (
33  const Specie& sp,
34  const table2D& rho
35 )
36 :
37  Specie(sp),
38  rho_(rho)
39 {}
40 
41 
42 template<class Specie>
44 (
45  const word& name,
46  const rhoTabulated<Specie>& ip
47 )
48 :
49  Specie(name, ip),
50  rho_(ip.rho_)
51 {}
52 
53 
54 template<class Specie>
57 {
59  (
60  new rhoTabulated<Specie>(*this)
61  );
62 }
63 
64 
65 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
66 
67 template<class Specie>
69 (
70  scalar p,
71  scalar T
72 ) const
73 {
74  return rho_.value(p, T);
75 }
76 
77 
78 template<class Specie>
79 inline Foam::scalar Foam::rhoTabulated<Specie>::h
80 (
81  scalar p,
82  scalar T
83 ) const
84 {
86  return 0;
87 }
88 
89 
90 template<class Specie>
91 inline Foam::scalar Foam::rhoTabulated<Specie>::Cp
92 (
93  scalar p,
94  scalar T
95 ) const
96 {
98  return 0;
99 }
100 
101 
102 template<class Specie>
103 inline Foam::scalar Foam::rhoTabulated<Specie>::e
104 (
105  scalar p,
106  scalar T
107 ) const
108 {
110  return 0;
111 }
112 
113 
114 template<class Specie>
116 (
117  scalar p,
118  scalar T
119 ) const
120 {
122  return 0;
123 }
124 
125 
126 template<class Specie>
128 (
129  scalar p,
130  scalar T
131 ) const
132 {
134  return 0;
135 }
136 
137 
138 template<class Specie>
140 (
141  scalar p,
142  scalar T
143 ) const
144 {
146  return 0;
147 }
148 
149 
150 template<class Specie>
152 (
153  scalar p,
154  scalar T
155 ) const
156 {
157  return rho_.dfdp(p, T);
158 }
159 
160 
161 template<class Specie>
162 inline Foam::scalar Foam::rhoTabulated<Specie>::Z
163 (
164  scalar p,
165  scalar T
166 ) const
167 {
168  return p/(rho(p, T)*this->R()*T);
169 }
170 
171 
172 template<class Specie>
174 (
175  scalar p,
176  scalar T
177 ) const
178 {
180  return 0;
181 }
182 
183 
184 template<class Specie>
186 (
187  scalar p,
188  scalar T
189 ) const
190 {
191  return -rho_.dfdT(p, T)/rho(p, T);
192 }
193 
194 
195 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
196 
197 template<class Specie>
199 (
200  const rhoTabulated<Specie>& pf
201 )
202 {
204 }
205 
206 
207 template<class Specie>
208 inline void Foam::rhoTabulated<Specie>::operator*=(const scalar s)
209 {
211 }
212 
213 
214 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
215 
216 template<class Specie>
217 inline Foam::rhoTabulated<Specie> Foam::operator+
218 (
219  const rhoTabulated<Specie>& pf1,
220  const rhoTabulated<Specie>& pf2
221 )
222 {
224  return pf1;
225 }
226 
227 
228 template<class Specie>
229 inline Foam::rhoTabulated<Specie> Foam::operator*
230 (
231  const scalar s,
232  const rhoTabulated<Specie>& pf
233 )
234 {
235  noCoefficientMixing(rhoTabulated);
236  return pf;
237 }
238 
239 
240 template<class Specie>
241 inline Foam::rhoTabulated<Specie> Foam::operator==
242 (
243  const rhoTabulated<Specie>& pf1,
244  const rhoTabulated<Specie>& pf2
245 )
246 {
247  noCoefficientMixing(rhoTabulated);
248  return pf1;
249 }
250 
251 
252 // ************************************************************************* //
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Incompressible of equation of state using uniform tabulated density vs pressure and temperature.
Definition: rhoTabulated.H:135
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].
rhoTabulated(const Specie &sp, const table2D &rho)
Construct from components.
Definition: rhoTabulatedI.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].
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
Definition: rhoTabulatedI.H:69
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: rhoTabulatedI.H:80
autoPtr< rhoTabulated > clone() const
Construct and return a clone.
Definition: rhoTabulatedI.H:56
scalar Cp(scalar p, scalar T) const
Return Cp contribution [J/(kg K].
Definition: rhoTabulatedI.H:92
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 [].
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))
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
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
#define noCoefficientMixing(Type)
Definition: specie.H:159