FreeFOAM The Cross-Platform CFD Toolkit
liquidMixture.C
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) 1991-2010 OpenCFD Ltd.
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 "liquidMixture.H"
27 #include <OpenFOAM/dictionary.H>
28 #include <specie/specie.H>
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 const Foam::scalar Foam::liquidMixture::TrMax = 0.999;
33 
34 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35 
37 (
38  const dictionary& thermophysicalProperties
39 )
40 :
41  components_(thermophysicalProperties.lookup("liquidComponents")),
42  properties_(components_.size())
43 {
44  // use sub-dictionary "liquidProperties" if possible to avoid
45  // collisions with identically named gas-phase entries
46  // (eg, H2O liquid vs. gas)
47  forAll(components_, i)
48  {
49  const dictionary* subDictPtr = thermophysicalProperties.subDictPtr
50  (
51  "liquidProperties"
52  );
53 
54  if (subDictPtr)
55  {
56  properties_.set
57  (
58  i,
59  liquid::New(subDictPtr->lookup(components_[i]))
60  );
61  }
62  else
63  {
64  properties_.set
65  (
66  i,
67  liquid::New(thermophysicalProperties.lookup(components_[i]))
68  );
69  }
70  }
71 }
72 
73 
74 // * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
75 
77 (
78  const dictionary& thermophysicalProperties
79 )
80 {
81  return autoPtr<liquidMixture>(new liquidMixture(thermophysicalProperties));
82 }
83 
84 
85 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
86 
87 Foam::scalar Foam::liquidMixture::Tc
88 (
89  const scalarField& x
90 ) const
91 {
92 
93  scalar vTc = 0.0;
94  scalar vc = 0.0;
95 
96  forAll(properties_, i)
97  {
98  scalar x1 = x[i]*properties_[i].Vc();
99  vc += x1;
100  vTc += x1*properties_[i].Tc();
101  }
102 
103  return vTc/vc;
104 }
105 
106 
107 Foam::scalar Foam::liquidMixture::Tpc
108 (
109  const scalarField& x
110 ) const
111 {
112  scalar Tpc = 0.0;
113  forAll(properties_, i)
114  {
115  Tpc += x[i]*properties_[i].Tc();
116  }
117 
118  return Tpc;
119 }
120 
121 
122 Foam::scalar Foam::liquidMixture::Ppc
123 (
124  const scalarField& x
125 ) const
126 {
127  scalar Vc = 0.0;
128  scalar Zc = 0.0;
129  forAll(properties_, i)
130  {
131  Vc += x[i]*properties_[i].Vc();
132  Zc += x[i]*properties_[i].Zc();
133  }
134 
135  return specie::RR*Zc*Tpc(x)/Vc;
136 }
137 
138 
139 Foam::scalar Foam::liquidMixture::omega
140 (
141  const scalarField& x
142 ) const
143 {
144  scalar omega = 0.0;
145  forAll(properties_, i)
146  {
147  omega += x[i]*properties_[i].omega();
148  }
149 
150  return omega;
151 }
152 
153 
155 (
156  const scalar p,
157  const scalar Tg,
158  const scalar Tl,
159  const scalarField& xg,
160  const scalarField& xl
161 ) const
162 {
163  scalarField xs(xl.size(), 0.0);
164 
165  // Raoult's Law
166  forAll(xs, i)
167  {
168  scalar Ti = min(TrMax*properties_[i].Tc(), Tl);
169  xs[i] = properties_[i].pv(p, Ti)*xl[i]/p;
170  }
171  return xs;
172 }
173 
174 
175 Foam::scalar Foam::liquidMixture::W
176 (
177  const scalarField& x
178 ) const
179 {
180  scalar W = 0.0;
181  forAll(properties_, i)
182  {
183  W += x[i]*properties_[i].W();
184  }
185 
186  return W;
187 }
188 
189 
191 (
192  const scalarField& X
193 ) const
194 {
195  scalarField Y(X/W(X));
196 
197  forAll(Y, i)
198  {
199  Y[i] *= properties_[i].W();
200  }
201 
202  return Y;
203 }
204 
205 
207 (
208  const scalarField& Y
209 ) const
210 {
211  scalarField X(Y.size());
212  scalar Winv = 0.0;
213  forAll(X, i)
214  {
215  Winv += Y[i]/properties_[i].W();
216  X[i] = Y[i]/properties_[i].W();
217  }
218 
219  return X/Winv;
220 }
221 
222 
223 Foam::scalar Foam::liquidMixture::rho
224 (
225  const scalar p,
226  const scalar T,
227  const scalarField& x
228 ) const
229 {
230  scalar v = 0.0;
231 
232  forAll(properties_, i)
233  {
234  if (x[i] > SMALL)
235  {
236  scalar Ti = min(TrMax*properties_[i].Tc(), T);
237  scalar rho = SMALL + properties_[i].rho(p, Ti);
238  v += x[i]*properties_[i].W()/rho;
239  }
240  }
241 
242  return W(x)/v;
243 }
244 
245 
246 Foam::scalar Foam::liquidMixture::pv
247 (
248  const scalar p,
249  const scalar T,
250  const scalarField& x
251 ) const
252 {
253  scalar pv = 0.0;
254 
255  forAll(properties_, i)
256  {
257  if (x[i] > SMALL)
258  {
259  scalar Ti = min(TrMax*properties_[i].Tc(), T);
260  pv += x[i]*properties_[i].pv(p, Ti)*properties_[i].W();
261  }
262  }
263 
264  return pv/W(x);
265 }
266 
267 
268 Foam::scalar Foam::liquidMixture::hl
269 (
270  const scalar p,
271  const scalar T,
272  const scalarField& x
273 ) const
274 {
275  scalar hl = 0.0;
276 
277  forAll(properties_, i)
278  {
279  if (x[i] > SMALL)
280  {
281  scalar Ti = min(TrMax*properties_[i].Tc(), T);
282  hl += x[i]*properties_[i].hl(p, Ti)*properties_[i].W();
283  }
284  }
285 
286  return hl/W(x);
287 }
288 
289 
290 Foam::scalar Foam::liquidMixture::cp
291 (
292  const scalar p,
293  const scalar T,
294  const scalarField& x
295 ) const
296 {
297  scalar cp = 0.0;
298 
299  forAll(properties_, i)
300  {
301  if (x[i] > SMALL)
302  {
303  scalar Ti = min(TrMax*properties_[i].Tc(), T);
304  cp += x[i]*properties_[i].cp(p, Ti)*properties_[i].W();
305  }
306  }
307 
308  return cp/W(x);
309 }
310 
311 
312 Foam::scalar Foam::liquidMixture::sigma
313 (
314  const scalar p,
315  const scalar T,
316  const scalarField& x
317 ) const
318 {
319  // sigma is based on surface mole fractions
320  // which is estimated from Raoult's Law
321  scalar sigma = 0.0;
322  scalarField Xs(x.size(), 0.0);
323  scalar XsSum = 0.0;
324  forAll(properties_, i)
325  {
326  scalar Ti = min(TrMax*properties_[i].Tc(), T);
327  scalar Pvs = properties_[i].pv(p, Ti);
328  scalar xs = x[i]*Pvs/p;
329  XsSum += xs;
330  Xs[i] = xs;
331  }
332 
333  forAll(properties_, i)
334  {
335  if (Xs[i] > SMALL)
336  {
337  scalar Ti = min(TrMax*properties_[i].Tc(), T);
338  sigma += (Xs[i]/XsSum)*properties_[i].sigma(p, Ti);
339  }
340  }
341 
342  return sigma;
343 }
344 
345 
346 Foam::scalar Foam::liquidMixture::mu
347 (
348  const scalar p,
349  const scalar T,
350  const scalarField& x
351 ) const
352 {
353  scalar mu = 0.0;
354 
355  forAll(properties_, i)
356  {
357  if (x[i] > SMALL)
358  {
359  scalar Ti = min(TrMax*properties_[i].Tc(), T);
360  mu += x[i]*log(properties_[i].mu(p, Ti));
361  }
362  }
363 
364  return exp(mu);
365 }
366 
367 
368 Foam::scalar Foam::liquidMixture::K
369 (
370  const scalar p,
371  const scalar T,
372  const scalarField& x
373 ) const
374 {
375  // calculate superficial volume fractions phii
376  scalarField phii(x.size(), 0.0);
377  scalar pSum = 0.0;
378 
379  forAll(properties_, i)
380  {
381  scalar Ti = min(TrMax*properties_[i].Tc(), T);
382 
383  scalar Vi = properties_[i].W()/properties_[i].rho(p, Ti);
384  phii[i] = x[i]*Vi;
385  pSum += phii[i];
386  }
387 
388  forAll(phii, i)
389  {
390  phii[i] /= pSum;
391  }
392 
393  scalar K = 0.0;
394 
395  forAll(properties_, i)
396  {
397  scalar Ti = min(TrMax*properties_[i].Tc(), T);
398 
399  forAll(properties_, j)
400  {
401  scalar Tj = min(TrMax*properties_[j].Tc(), T);
402 
403  scalar Kij =
404  2.0
405  /(
406  1.0/properties_[i].K(p, Ti)
407  + 1.0/properties_[j].K(p, Tj)
408  );
409  K += phii[i]*phii[j]*Kij;
410  }
411  }
412 
413  return K;
414 }
415 
416 
417 Foam::scalar Foam::liquidMixture::D
418 (
419  const scalar p,
420  const scalar T,
421  const scalarField& x
422 ) const
423 {
424  // Blanc's law
425  scalar Dinv = 0.0;
426 
427  forAll(properties_, i)
428  {
429  if (x[i] > SMALL)
430  {
431  scalar Ti = min(TrMax*properties_[i].Tc(), T);
432  Dinv += x[i]/properties_[i].D(p, Ti);
433  }
434  }
435 
436  return 1.0/Dinv;
437 }
438 
439 
440 // ************************ vim: set sw=4 sts=4 et: ************************ //