rhoConstI.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) 2012-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 "rhoConst.H"
29
30// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
31
32template<class Specie>
34(
35 const Specie& sp,
36 const scalar rho
37)
38:
39 Specie(sp),
40 rho_(rho)
41{}
42
43
44// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45
46template<class Specie>
48(
49 const word& name,
50 const rhoConst<Specie>& ico
51)
52:
53 Specie(name, ico),
54 rho_(ico.rho_)
55{}
56
57
58template<class Specie>
61{
62 return autoPtr<rhoConst<Specie>>::New(*this);
63}
64
65
66// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
67
68template<class Specie>
69inline Foam::scalar Foam::rhoConst<Specie>::rho(scalar p, scalar T) const
70{
71 return rho_;
72}
73
74
75template<class Specie>
76inline Foam::scalar Foam::rhoConst<Specie>::H(scalar p, scalar T) const
77{
78 return 0;
79}
80
81
82template<class Specie>
83inline Foam::scalar Foam::rhoConst<Specie>::Cp(scalar p, scalar T) const
84{
85 return 0;
86}
87
88
89template<class Specie>
90inline Foam::scalar Foam::rhoConst<Specie>::E(scalar p, scalar T) const
91{
92 return 0;
93}
94
95
96template<class Specie>
97inline Foam::scalar Foam::rhoConst<Specie>::Cv(scalar p, scalar T) const
98{
99 return 0;
100}
101
102
103template<class Specie>
104inline Foam::scalar Foam::rhoConst<Specie>::S(scalar p, scalar T) const
105{
106 return 0;
107}
108
109
110template<class Specie>
111inline Foam::scalar Foam::rhoConst<Specie>::psi(scalar p, scalar T) const
112{
113 return 0;
114}
115
116
117template<class Specie>
118inline Foam::scalar Foam::rhoConst<Specie>::Z(scalar p, scalar T) const
119{
120 return 0;
121}
122
123
124template<class Specie>
125inline Foam::scalar Foam::rhoConst<Specie>::CpMCv(scalar p, scalar T) const
126{
127 return 0;
128}
129
130
131// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
132
133template<class Specie>
135{
136 scalar Y1 = this->Y();
137 Specie::operator+=(ico);
138
139 if (mag(this->Y()) > SMALL)
140 {
141 Y1 /= this->Y();
142 const scalar Y2 = ico.Y()/this->Y();
143
144 rho_ = Y1*rho_ + Y2*ico.rho_;
145 }
146}
147
148
149template<class Specie>
150inline void Foam::rhoConst<Specie>::operator*=(const scalar s)
151{
152 Specie::operator*=(s);
153}
154
155
156// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
157
158template<class Specie>
159inline Foam::rhoConst<Specie> Foam::operator+
160(
161 const rhoConst<Specie>& ico1,
162 const rhoConst<Specie>& ico2
163)
164{
165 Specie sp
166 (
167 static_cast<const Specie&>(ico1)
168 + static_cast<const Specie&>(ico2)
169 );
170
171 if (mag(sp.Y()) < SMALL)
172 {
173 return rhoConst<Specie>
174 (
175 sp,
176 ico1.rho_
177 );
178 }
179 else
180 {
181 const scalar Y1 = ico1.Y()/sp.Y();
182 const scalar Y2 = ico2.Y()/sp.Y();
183
184 return rhoConst<Specie>
185 (
186 sp,
187 Y1*ico1.rho_ + Y2*ico2.rho_
188 );
189 }
190}
191
192
193template<class Specie>
194inline Foam::rhoConst<Specie> Foam::operator*
195(
196 const scalar s,
197 const rhoConst<Specie>& ico
198)
199{
200 return rhoConst<Specie>(s*static_cast<const Specie&>(ico), ico.rho_);
201}
202
203
204template<class Specie>
205inline Foam::rhoConst<Specie> Foam::operator==
206(
207 const rhoConst<Specie>& ico1,
208 const rhoConst<Specie>& ico2
209)
210{
211 Specie sp
212 (
213 static_cast<const Specie&>(ico1)
214 == static_cast<const Specie&>(ico2)
215 );
216
217 const scalar Y1 = ico1.Y()/sp.Y();
218 const scalar Y2 = ico2.Y()/sp.Y();
219
220 return rhoConst<Specie>
221 (
222 sp,
223 Y2*ico2.rho_ - Y1*ico1.rho_
224 );
225}
226
227
228// ************************************************************************* //
compactSpatialTensor S
The joint motion sub-space (3-DoF)
Definition: joint.H:132
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
tmp< GeometricField< Type, faPatchField, areaMesh > > H() const
Return the H operation source.
Definition: faMatrix.C:633
RhoConst (rho = const) of state.
Definition: rhoConst.H:91
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
Definition: rhoConstI.H:125
autoPtr< rhoConst > clone() const
Construct and return a clone.
Definition: rhoConstI.H:60
void operator*=(const scalar)
Definition: rhoConstI.H:150
void operator+=(const rhoConst &)
Definition: rhoConstI.H:134
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
PtrList< volScalarField > & Y
const volScalarField & psi
const volScalarField & T
const volScalarField & Cv
Definition: EEqn.H:8
const volScalarField & Cp
Definition: EEqn.H:7
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))
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59