rPolynomialI.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) 2019-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 "rPolynomial.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class Specie>
32 (
33  const Specie& sp,
34  const coeffList& coeffs
35 )
36 :
37  Specie(sp),
38  C_(coeffs)
39 {}
40 
41 
42 template<class Specie>
44 (
45  const word& name,
46  const rPolynomial<Specie>& rp
47 )
48 :
49  Specie(name, rp),
50  C_(rp.C_)
51 {}
52 
53 
54 template<class Specie>
57 {
59 }
60 
61 
62 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
63 
64 template<class Specie>
65 inline Foam::scalar Foam::rPolynomial<Specie>::rho(scalar p, scalar T) const
66 {
67  return 1/(C_[0] + (C_[1] + C_[2]*T - C_[4]*p)*T - C_[3]*p);
68 }
69 
70 
71 template<class Specie>
72 inline Foam::scalar Foam::rPolynomial<Specie>::h(scalar p, scalar T) const
73 {
74  return 0;
75 }
76 
77 
78 template<class Specie>
79 inline Foam::scalar Foam::rPolynomial<Specie>::Cp(scalar p, scalar T) const
80 {
81  return 0;
82 }
83 
84 
85 template<class Specie>
86 inline Foam::scalar Foam::rPolynomial<Specie>::e(scalar p, scalar T) const
87 {
88  return 0;
89 }
90 
91 
92 template<class Specie>
93 inline Foam::scalar Foam::rPolynomial<Specie>::Cv(scalar p, scalar T) const
94 {
95  return 0;
96 }
97 
98 
99 template<class Specie>
100 inline Foam::scalar Foam::rPolynomial<Specie>::sp(scalar p, scalar T) const
101 {
102  return 0;
103 }
104 
105 
106 template<class Specie>
107 inline Foam::scalar Foam::rPolynomial<Specie>::sv(scalar p, scalar T) const
108 {
110  return 0;
111 }
112 
113 
114 template<class Specie>
115 inline Foam::scalar Foam::rPolynomial<Specie>::psi(scalar p, scalar T) const
116 {
117  return sqr(rho(p, T))*(C_[3] + C_[4]*T);
118 }
119 
120 
121 template<class Specie>
122 inline Foam::scalar Foam::rPolynomial<Specie>::Z(scalar p, scalar T) const
123 {
124  return p/(rho(p, T)*this->R()*T);
125 }
126 
127 
128 template<class Specie>
129 inline Foam::scalar Foam::rPolynomial<Specie>::CpMCv(scalar p, scalar T) const
130 {
131  return 0;
132 }
133 
134 
135 template<class Specie>
136 inline Foam::scalar Foam::rPolynomial<Specie>::alphav(scalar p, scalar T) const
137 {
138  return this->rho(p, T)*(C_[1] + 2*C_[2]*T - C_[4]*p);
139 }
140 
141 
142 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
143 
144 template<class Specie>
146 (
147  const rPolynomial<Specie>& rp
148 )
149 {
150  const scalar Y1 = this->Y();
151  Specie::operator+=(rp);
152 
153  if (mag(this->Y()) > small)
154  {
155  C_ = (Y1*C_ + rp.Y()*rp.C_)/this->Y();
156  }
157 }
158 
159 
160 template<class Specie>
161 inline void Foam::rPolynomial<Specie>::operator*=(const scalar s)
162 {
163  Specie::operator*=(s);
164 }
165 
166 
167 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
168 
169 template<class Specie>
170 inline Foam::rPolynomial<Specie> Foam::operator+
171 (
172  const rPolynomial<Specie>& rp1,
173  const rPolynomial<Specie>& rp2
174 )
175 {
176  Specie sp
177  (
178  static_cast<const Specie&>(rp1)
179  + static_cast<const Specie&>(rp2)
180  );
181 
182  if (mag(sp.Y()) < small)
183  {
184  return rPolynomial<Specie>
185  (
186  sp,
187  rp1.C_
188  );
189  }
190  else
191  {
192  return rPolynomial<Specie>
193  (
194  sp,
195  (rp1.Y()*rp1.C_ + rp2.Y()*rp2.C_)/sp.Y()
196  );
197  }
198 
199  return rp1;
200 }
201 
202 
203 template<class Specie>
204 inline Foam::rPolynomial<Specie> Foam::operator*
205 (
206  const scalar s,
207  const rPolynomial<Specie>& rp
208 )
209 {
210  return rPolynomial<Specie>
211  (
212  s*static_cast<const Specie&>(rp),
213  rp.C_
214  );
215 }
216 
217 
218 template<class Specie>
219 inline Foam::rPolynomial<Specie> Foam::operator==
220 (
221  const rPolynomial<Specie>& rp1,
222  const rPolynomial<Specie>& rp2
223 )
224 {
225  return rPolynomial<Specie>
226  (
227  static_cast<const Specie&>(rp1) == static_cast<const Specie&>(rp2),
228  rPolynomial<Specie>::coeffList::uniform(NaN)
229  );
230 }
231 
232 
233 // ************************************************************************* //
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Reciprocal polynomial equation of state for liquids and solids.
Definition: rPolynomial.H:118
scalar Cv(scalar p, scalar T) const
Return Cv contribution [J/(kg K].
Definition: rPolynomialI.H:93
scalar psi(scalar p, scalar T) const
Return compressibility [s^2/m^2].
Definition: rPolynomialI.H:115
scalar alphav(const scalar p, const scalar T) const
Return volumetric coefficient of thermal expansion [1/T].
Definition: rPolynomialI.H:136
scalar e(const scalar p, const scalar T) const
Return internal energy contribution [J/kg].
Definition: rPolynomialI.H:86
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
Definition: rPolynomialI.H:65
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
Definition: rPolynomialI.H:129
scalar h(const scalar p, const scalar T) const
Return enthalpy contribution [J/kg].
Definition: rPolynomialI.H:72
scalar Cp(scalar p, scalar T) const
Return Cp contribution [J/(kg K].
Definition: rPolynomialI.H:79
scalar sv(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cv/T [J/kg/K].
Definition: rPolynomialI.H:107
rPolynomial(const Specie &sp, const coeffList &coeffs)
Construct from components.
Definition: rPolynomialI.H:32
scalar sp(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cp/T [J/kg/K].
Definition: rPolynomialI.H:100
scalar Z(scalar p, scalar T) const
Return compression factor [].
Definition: rPolynomialI.H:122
autoPtr< rPolynomial > clone() const
Construct and return a clone.
Definition: rPolynomialI.H:56
void operator*=(const scalar)
Definition: rPolynomialI.H:161
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))
dimensionedSymmTensor sqr(const dimensionedVector &dv)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
dimensioned< scalar > mag(const dimensioned< Type > &)
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
PtrList< volScalarField > & Y