coefficientWilkeMultiComponentMixture.C
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) 2020 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 
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template<class ThermoType>
33 (
34  const dictionary& thermoDict,
35  const fvMesh& mesh,
36  const word& phaseName
37 )
38 :
40  (
41  thermoDict,
42  mesh,
43  phaseName
44  ),
45  mixture_("mixture", this->specieThermos()[0]),
46  transportMixture_(this->specieThermos())
47 {}
48 
49 
50 template<class ThermoType>
53 (
54  const PtrList<ThermoType>& specieThermos
55 )
56 :
57  specieThermos_(specieThermos),
58  M_(specieThermos.size()),
59  A_(specieThermos.size()),
60  B_(specieThermos.size()),
61  X_(specieThermos.size()),
62  mu_(specieThermos.size()),
63  w_(specieThermos.size()),
64  muCached_(false)
65 {
66  forAll(specieThermos_, i)
67  {
68  M_[i] = specieThermos[i].W();
69  }
70 
71  forAll(M_, i)
72  {
73  forAll(M_, j)
74  {
75  if (i != j)
76  {
77  A_(i, j) = ((4/sqrt(2.0))*sqrt(1 + M_[i]/M_[j]));
78  B_(i, j) = sqrt(M_[j]/M_[i]);
79  }
80  }
81  }
82 }
83 
84 
85 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
86 
87 template<class ThermoType>
90 (
91  scalar p,
92  scalar T
93 ) const
94 {
95  forAll(mu_, i)
96  {
97  mu_[i] = specieThermos_[i].mu(p, T);
98  }
99 
100  forAll(M_, i)
101  {
102  scalar sumXphi = 0;
103 
104  forAll(M_, j)
105  {
106  if (i != j)
107  {
108  const scalar phiij =
109  sqr(1 + sqrt((mu_[i]/mu_[j])*B_(i, j)))/A_(i, j);
110 
111  sumXphi += X_[j]*phiij;
112  }
113  else
114  {
115  sumXphi += X_[j];
116  }
117  }
118 
119  w_[i] = X_[i]/sumXphi;
120  }
121 }
122 
123 
124 template<class ThermoType>
125 Foam::scalar
127 mu
128 (
129  scalar p,
130  scalar T
131 ) const
132 {
133  WilkeWeights(p, T);
134 
135  scalar mu = 0;
136  forAll(w_, i)
137  {
138  mu += w_[i]*mu_[i];
139  }
140 
141  return mu;
142 }
143 
144 
145 template<class ThermoType>
146 Foam::scalar
148 kappa
149 (
150  scalar p,
151  scalar T
152 ) const
153 {
154  if (!muCached_)
155  {
156  WilkeWeights(p, T);
157  }
158 
159  scalar kappa = 0;
160  forAll(w_, i)
161  {
162  kappa += w_[i]*specieThermos_[i].kappa(p, T);
163  }
164 
165  return kappa;
166 }
167 
168 
169 template<class ThermoType>
170 const typename
173 (
174  const label celli
175 ) const
176 {
177  mixture_ = this->Y()[0][celli]*this->specieThermos()[0];
178 
179  for (label i=1; i<this->Y().size(); i++)
180  {
181  mixture_ += this->Y()[i][celli]*this->specieThermos()[i];
182  }
183 
184  return mixture_;
185 }
186 
187 
188 template<class ThermoType>
189 const typename
192 (
193  const label patchi,
194  const label facei
195 ) const
196 {
197  mixture_ =
198  this->Y()[0].boundaryField()[patchi][facei]
199  *this->specieThermos()[0];
200 
201  for (label i=1; i<this->Y().size(); i++)
202  {
203  mixture_ +=
204  this->Y()[i].boundaryField()[patchi][facei]
205  *this->specieThermos()[i];
206  }
207 
208  return mixture_;
209 }
210 
211 
212 template<class ThermoType>
213 const typename
216 (
217  const label celli
218 ) const
219 {
220  transportMixture_.muCached_ = false;
221 
222  scalarList& X = transportMixture_.X_;
223 
224  scalar sumX = 0;
225 
226  forAll(X, i)
227  {
228  X[i] = this->Y()[i][celli]/this->specieThermos()[i].W();
229  sumX += X[i];
230  }
231 
232  forAll(X, i)
233  {
234  X[i] /= sumX;
235  }
236 
237  return transportMixture_;
238 }
239 
240 
241 template<class ThermoType>
242 const typename
246 (
247  const label patchi,
248  const label facei
249 ) const
250 {
251  transportMixture_.muCached_ = false;
252 
253  scalarList& X = transportMixture_.X_;
254 
255  scalar sumX = 0;
256 
257  forAll(X, i)
258  {
259  X[i] =
260  this->Y()[i].boundaryField()[patchi][facei]
261  /this->specieThermos()[i].W();
262  sumX += X[i];
263  }
264 
265  forAll(X, i)
266  {
267  X[i] /= sumX;
268  }
269 
270  return transportMixture_;
271 }
272 
273 template<class ThermoType>
274 const typename
278 (
279  const label celli,
280  const thermoMixtureType& thermoMixture
281 ) const
282 {
283  cellTransportMixture(celli);
284  transportMixture_.muCached_ = true;
285  return transportMixture_;
286 }
287 
288 
289 template<class ThermoType>
290 const typename
294 (
295  const label patchi,
296  const label facei,
297  const thermoMixtureType& thermoMixture
298 ) const
299 {
300  patchFaceTransportMixture(patchi, facei);
301  transportMixture_.muCached_ = true;
302  return transportMixture_;
303 }
304 
305 
306 // ************************************************************************* //
const thermoMixtureType & cellThermoMixture(const label celli) const
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Thermophysical properties mixing class which applies mass-fraction weighted mixing to the thermodynam...
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const thermoMixtureType & patchFaceThermoMixture(const label patchi, const label facei) const
dimensionedScalar sqrt(const dimensionedScalar &ds)
fvMesh & mesh
const transportMixtureType & cellTransportMixture(const label celli) const
const PtrList< ThermoType > & specieThermos() const
Return the raw specie thermodynamic data.
A class for handling words, derived from string.
Definition: word.H:59
coefficientWilkeMultiComponentMixture(const dictionary &, const fvMesh &, const word &)
Construct from dictionary, mesh and phase name.
Foam::multiComponentMixture.
scalar mu(const scalar p, const scalar T) const
Dynamic viscosity [kg/m/s].
const transportMixtureType & patchFaceTransportMixture(const label patchi, const label facei) const
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label patchi
scalar kappa(const scalar p, const scalar T) const
Thermal conductivity [W/m/K].
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
ThermoType::thermoType thermoMixtureType
Mixing type for thermodynamic properties.