constTransportI.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 | Copyright (C) 2011-2017 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 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
27 
28 template<class Thermo>
30 (
31  const Thermo& t,
32  const scalar mu,
33  const scalar Pr
34 )
35 :
36  Thermo(t),
37  mu_(mu),
38  rPr_(1.0/Pr)
39 {}
40 
41 
42 template<class Thermo>
44 (
45  const word& name,
46  const constTransport& ct
47 )
48 :
49  Thermo(name, ct),
50  mu_(ct.mu_),
51  rPr_(ct.rPr_)
52 {}
53 
54 
55 template<class Thermo>
58 {
60  (
61  new constTransport<Thermo>(*this)
62  );
63 }
64 
65 
66 template<class Thermo>
69 (
70  const dictionary& dict
71 )
72 {
74  (
76  );
77 }
78 
79 
80 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
81 
82 template<class Thermo>
83 inline Foam::scalar Foam::constTransport<Thermo>::mu
84 (
85  const scalar p,
86  const scalar T
87 ) const
88 {
89  return mu_;
90 }
91 
92 
93 template<class Thermo>
94 inline Foam::scalar Foam::constTransport<Thermo>::kappa
95 (
96  const scalar p,
97  const scalar T
98 ) const
99 {
100  return this->Cp(p, T)*mu(p, T)*rPr_;
101 }
102 
103 
104 template<class Thermo>
105 inline Foam::scalar Foam::constTransport<Thermo>::alphah
106 (
107  const scalar p,
108  const scalar T
109 ) const
110 {
111  return mu(p, T)*rPr_;
112 }
113 
114 
115 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
116 
117 template<class Thermo>
118 inline void Foam::constTransport<Thermo>::operator=
119 (
120  const constTransport<Thermo>& ct
121 )
122 {
123  Thermo::operator=(ct);
124 
125  mu_ = ct.mu_;
126  rPr_ = ct.rPr_;
127 }
128 
129 
130 template<class Thermo>
131 inline void Foam::constTransport<Thermo>::operator+=
132 (
133  const constTransport<Thermo>& st
134 )
135 {
136  scalar Y1 = this->Y();
137 
138  Thermo::operator+=(st);
139 
140  if (mag(this->Y()) > SMALL)
141  {
142  Y1 /= this->Y();
143  scalar Y2 = st.Y()/this->Y();
144 
145  mu_ = Y1*mu_ + Y2*st.mu_;
146  rPr_ = 1.0/(Y1/rPr_ + Y2/st.rPr_);
147  }
148 }
149 
150 
151 template<class Thermo>
152 inline void Foam::constTransport<Thermo>::operator*=
153 (
154  const scalar s
155 )
156 {
157  Thermo::operator*=(s);
158 }
159 
160 
161 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
162 
163 template<class Thermo>
164 inline Foam::constTransport<Thermo> Foam::operator+
165 (
166  const constTransport<Thermo>& ct1,
167  const constTransport<Thermo>& ct2
168 )
169 {
170  Thermo t
171  (
172  static_cast<const Thermo&>(ct1) + static_cast<const Thermo&>(ct2)
173  );
174 
175  if (mag(t.Y()) < SMALL)
176  {
178  (
179  t,
180  0,
181  ct1.rPr_
182  );
183  }
184  else
185  {
186  scalar Y1 = ct1.Y()/t.Y();
187  scalar Y2 = ct2.Y()/t.Y();
188 
190  (
191  t,
192  Y1*ct1.mu_ + Y2*ct2.mu_,
193  1.0/(Y1/ct1.rPr_ + Y2/ct2.rPr_)
194  );
195  }
196 }
197 
198 
199 template<class Thermo>
200 inline Foam::constTransport<Thermo> Foam::operator*
201 (
202  const scalar s,
203  const constTransport<Thermo>& ct
204 )
205 {
207  (
208  s*static_cast<const Thermo&>(ct),
209  ct.mu_,
210  1.0/ct.rPr_
211  );
212 }
213 
214 
215 // ************************************************************************* //
PtrList< volScalarField > & Y1
Definition: YEqns.H:8
dimensionedScalar Pr("Pr", dimless, laminarTransport)
dictionary dict
scalar kappa(const scalar p, const scalar T) const
Thermal conductivity [W/mK].
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Constant properties Transport package. Templated into a given thermodynamics package (needed for ther...
PtrList< volScalarField > & Y2
Definition: YEqns.H:9
scalar alphah(const scalar p, const scalar T) const
Thermal diffusivity of enthalpy [kg/ms].
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.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
A class for handling words, derived from string.
Definition: word.H:59
autoPtr< constTransport > clone() const
Construct and return a clone.
const dimensionedScalar mu
Atomic mass unit.
PtrList< volScalarField > & Y
dimensioned< scalar > mag(const dimensioned< Type > &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
static autoPtr< constTransport > New(const dictionary &dict)
scalar mu(const scalar p, const scalar T) const
Dynamic viscosity [kg/ms].