sutherlandTransportI.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) 2011-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 "specie.H"
27 
28 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
29 
30 template<class Thermo>
32 (
33  const scalar mu1,
34  const scalar T1,
35  const scalar mu2,
36  const scalar T2
37 )
38 {
39  scalar rootT1 = sqrt(T1);
40  scalar mu1rootT2 = mu1*sqrt(T2);
41  scalar mu2rootT1 = mu2*rootT1;
42 
43  Ts_ = (mu2rootT1 - mu1rootT2)/(mu1rootT2/T1 - mu2rootT1/T2);
44 
45  As_ = mu1*(1.0 + Ts_/T1)/rootT1;
46 }
47 
48 
49 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
50 
51 template<class Thermo>
53 (
54  const Thermo& t,
55  const scalar As,
56  const scalar Ts
57 )
58 :
59  Thermo(t),
60  As_(As),
61  Ts_(Ts)
62 {}
63 
64 
65 template<class Thermo>
67 (
68  const Thermo& t,
69  const scalar mu1,
70  const scalar T1,
71  const scalar mu2,
72  const scalar T2
73 )
74 :
75  Thermo(t)
76 {
77  calcCoeffs(mu1, T1, mu2, T2);
78 }
79 
80 
81 template<class Thermo>
83 (
84  const word& name,
85  const sutherlandTransport& st
86 )
87 :
88  Thermo(name, st),
89  As_(st.As_),
90  Ts_(st.Ts_)
91 {}
92 
93 
94 template<class Thermo>
97 {
99  (
100  new sutherlandTransport<Thermo>(*this)
101  );
102 }
103 
104 
105 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
106 
107 template<class Thermo>
109 (
110  const scalar p,
111  const scalar T
112 ) const
113 {
114  return As_*::sqrt(T)/(1.0 + Ts_/T);
115 }
116 
117 
118 template<class Thermo>
120 (
121  const scalar p,
122  const scalar T
123 ) const
124 {
125  scalar Cv_ = this->Cv(p, T);
126  return mu(p, T)*Cv_*(1.32 + 1.77*this->R()/Cv_);
127 }
128 
129 
130 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
131 
132 template<class Thermo>
134 (
136 )
137 {
138  scalar Y1 = this->Y();
139 
140  Thermo::operator+=(st);
141 
142  if (mag(this->Y()) > small)
143  {
144  Y1 /= this->Y();
145  scalar Y2 = st.Y()/this->Y();
146 
147  As_ = Y1*As_ + Y2*st.As_;
148  Ts_ = Y1*Ts_ + Y2*st.Ts_;
149  }
150 }
151 
152 
153 template<class Thermo>
155 (
156  const scalar s
157 )
158 {
159  Thermo::operator*=(s);
160 }
161 
162 
163 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
164 
165 template<class Thermo>
166 inline Foam::sutherlandTransport<Thermo> Foam::operator+
167 (
168  const sutherlandTransport<Thermo>& st1,
169  const sutherlandTransport<Thermo>& st2
170 )
171 {
172  Thermo t
173  (
174  static_cast<const Thermo&>(st1) + static_cast<const Thermo&>(st2)
175  );
176 
177  if (mag(t.Y()) < small)
178  {
180  (
181  t,
182  0,
183  st1.As_,
184  st1.Ts_
185  );
186  }
187  else
188  {
189  scalar Y1 = st1.Y()/t.Y();
190  scalar Y2 = st2.Y()/t.Y();
191 
192  return sutherlandTransport<Thermo>
193  (
194  t,
195  Y1*st1.As_ + Y2*st2.As_,
196  Y1*st1.Ts_ + Y2*st2.Ts_
197  );
198  }
199 }
200 
201 
202 template<class Thermo>
203 inline Foam::sutherlandTransport<Thermo> Foam::operator*
204 (
205  const scalar s,
206  const sutherlandTransport<Thermo>& st
207 )
208 {
209  return sutherlandTransport<Thermo>
210  (
211  s*static_cast<const Thermo&>(st),
212  st.As_,
213  st.Ts_
214  );
215 }
216 
217 
218 // ************************************************************************* //
scalar Cv(const scalar p, const scalar T) const
Definition: HtoEthermo.H:2
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Transport package using Sutherland's formula for viscosity:
scalar mu(const scalar p, const scalar T) const
Dynamic viscosity [kg/m/s].
scalar kappa(const scalar p, const scalar T) const
Thermal conductivity [W/m/K].
sutherlandTransport(const Thermo &t, const scalar As, const scalar Ts)
Construct from components.
autoPtr< sutherlandTransport > clone() const
Construct and return a clone.
A class for handling words, derived from string.
Definition: word.H:62
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))
const dimensionedScalar mu
Atomic mass unit.
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
dimensionedScalar sqrt(const dimensionedScalar &ds)
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