sutherlandTransportI.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-2013 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 "specie.H"
27 
28 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
29 
30 template<class Thermo>
32 (
33  const scalar mu1, const scalar T1,
34  const scalar mu2, const scalar T2
35 )
36 {
37  scalar rootT1 = sqrt(T1);
38  scalar mu1rootT2 = mu1*sqrt(T2);
39  scalar mu2rootT1 = mu2*rootT1;
40 
41  Ts_ = (mu2rootT1 - mu1rootT2)/(mu1rootT2/T1 - mu2rootT1/T2);
42 
43  As_ = mu1*(1.0 + Ts_/T1)/rootT1;
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
49 template<class Thermo>
51 (
52  const Thermo& t,
53  const scalar As,
54  const scalar Ts
55 )
56 :
57  Thermo(t),
58  As_(As),
59  Ts_(Ts)
60 {}
61 
62 
63 template<class Thermo>
65 (
66  const Thermo& t,
67  const scalar mu1, const scalar T1,
68  const scalar mu2, const scalar T2
69 )
70 :
71  Thermo(t)
72 {
73  calcCoeffs(mu1, T1, mu2, T2);
74 }
75 
76 
77 template<class Thermo>
79 (
80  const word& name,
81  const sutherlandTransport& st
82 )
83 :
84  Thermo(name, st),
85  As_(st.As_),
86  Ts_(st.Ts_)
87 {}
88 
89 
90 template<class Thermo>
93 {
95  (
97  );
98 }
99 
100 
101 template<class Thermo>
104 (
105  Istream& is
106 )
107 {
109  (
111  );
112 }
113 
114 
115 template<class Thermo>
118 (
119  const dictionary& dict
120 )
121 {
123  (
125  );
126 }
127 
128 
129 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
130 
131 template<class Thermo>
132 inline Foam::scalar Foam::sutherlandTransport<Thermo>::mu
133 (
134  const scalar p,
135  const scalar T
136 ) const
137 {
138  return As_*::sqrt(T)/(1.0 + Ts_/T);
139 }
140 
141 
142 template<class Thermo>
144 (
145  const scalar p, const scalar T
146 ) const
147 {
148  scalar Cv_ = this->Cv(p, T);
149  return mu(p, T)*Cv_*(1.32 + 1.77*this->R()/Cv_);
150 }
151 
152 
153 template<class Thermo>
155 (
156  const scalar p,
157  const scalar T
158 ) const
159 {
160 
161  return kappa(p, T)/this->Cpv(p, T);
162 }
163 
164 
165 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
166 
167 template<class Thermo>
169 Foam::sutherlandTransport<Thermo>::operator=
170 (
172 )
173 {
174  Thermo::operator=(st);
175 
176  As_ = st.As_;
177  Ts_ = st.Ts_;
178 
179  return *this;
180 }
181 
182 
183 template<class Thermo>
184 inline void Foam::sutherlandTransport<Thermo>::operator+=
185 (
187 )
188 {
189  scalar molr1 = this->nMoles();
190 
191  Thermo::operator+=(st);
192 
193  molr1 /= this->nMoles();
194  scalar molr2 = st.nMoles()/this->nMoles();
195 
196  As_ = molr1*As_ + molr2*st.As_;
197  Ts_ = molr1*Ts_ + molr2*st.Ts_;
198 }
199 
200 
201 template<class Thermo>
202 inline void Foam::sutherlandTransport<Thermo>::operator-=
203 (
205 )
206 {
207  scalar molr1 = this->nMoles();
208 
209  Thermo::operator-=(st);
210 
211  molr1 /= this->nMoles();
212  scalar molr2 = st.nMoles()/this->nMoles();
213 
214  As_ = molr1*As_ - molr2*st.As_;
215  Ts_ = molr1*Ts_ - molr2*st.Ts_;
216 }
217 
218 
219 template<class Thermo>
220 inline void Foam::sutherlandTransport<Thermo>::operator*=
221 (
222  const scalar s
223 )
224 {
225  Thermo::operator*=(s);
226 }
227 
228 
229 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
230 
231 template<class Thermo>
232 inline Foam::sutherlandTransport<Thermo> Foam::operator+
233 (
234  const sutherlandTransport<Thermo>& st1,
235  const sutherlandTransport<Thermo>& st2
236 )
237 {
238  Thermo t
239  (
240  static_cast<const Thermo&>(st1) + static_cast<const Thermo&>(st2)
241  );
242 
243  scalar molr1 = st1.nMoles()/t.nMoles();
244  scalar molr2 = st2.nMoles()/t.nMoles();
245 
247  (
248  t,
249  molr1*st1.As_ + molr2*st2.As_,
250  molr1*st1.Ts_ + molr2*st2.Ts_
251  );
252 }
253 
254 
255 template<class Thermo>
256 inline Foam::sutherlandTransport<Thermo> Foam::operator-
257 (
258  const sutherlandTransport<Thermo>& st1,
259  const sutherlandTransport<Thermo>& st2
260 )
261 {
262  Thermo t
263  (
264  static_cast<const Thermo&>(st1) - static_cast<const Thermo&>(st2)
265  );
266 
267  scalar molr1 = st1.nMoles()/t.nMoles();
268  scalar molr2 = st2.nMoles()/t.nMoles();
269 
271  (
272  t,
273  molr1*st1.As_ - molr2*st2.As_,
274  molr1*st1.Ts_ - molr2*st2.Ts_
275  );
276 }
277 
278 
279 template<class Thermo>
280 inline Foam::sutherlandTransport<Thermo> Foam::operator*
281 (
282  const scalar s,
284 )
285 {
287  (
288  s*static_cast<const Thermo&>(st),
289  st.As_,
290  st.Ts_
291  );
292 }
293 
294 
295 template<class Thermo>
296 inline Foam::sutherlandTransport<Thermo> Foam::operator==
297 (
298  const sutherlandTransport<Thermo>& st1,
299  const sutherlandTransport<Thermo>& st2
300 )
301 {
302  return st2 - st1;
303 }
304 
305 
306 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
const dimensionedScalar mu
Atomic mass unit.
Definition: createFields.H:13
autoPtr< sutherlandTransport > clone() const
Construct and return a clone.
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 ))
#define R(A, B, C, D, E, F, K, M)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
A class for handling words, derived from string.
Definition: word.H:59
sutherlandTransport(const Thermo &t, const scalar As, const scalar Ts)
Construct from components.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const volScalarField & T
Definition: createFields.H:25
dictionary dict
static autoPtr< sutherlandTransport > New(Istream &is)
scalar alphah(const scalar p, const scalar T) const
Thermal diffusivity of enthalpy [kg/ms].
scalar kappa(const scalar p, const scalar T) const
Thermal conductivity [W/mK].
Transport package using Sutherland&#39;s formula.
scalar mu(const scalar p, const scalar T) const
Dynamic viscosity [kg/ms].
const dimensionedScalar kappa
Coulomb constant: default SI units: [N.m2/C2].
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117