OpenWAM
THTM_Fluids.h
1 /*--------------------------------------------------------------------------------*\
2 |==========================|
3  |\\ /\ /\ // O pen | OpenWAM: The Open Source 1D Gas-Dynamic Code
4  | \\ | X | // W ave |
5  | \\ \/_\/ // A ction | CMT-Motores Termicos / Universidad Politecnica Valencia
6  | \\/ \// M odel |
7  ----------------------------------------------------------------------------------
8  | License
9  |
10  | This file is part of OpenWAM.
11  |
12  | OpenWAM is free software: you can redistribute it and/or modify
13  | it under the terms of the GNU General Public License as published by
14  | the Free Software Foundation, either version 3 of the License, or
15  | (at your option) any later version.
16  |
17  | OpenWAM is distributed in the hope that it will be useful,
18  | but WITHOUT ANY WARRANTY; without even the implied warranty of
19  | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  | GNU General Public License for more details.
21  |
22  | You should have received a copy of the GNU General Public License
23  | along with OpenWAM. If not, see <http://www.gnu.org/licenses/>.
24  |
25  \*--------------------------------------------------------------------------------*/
26 
27 // ---------------------------------------------------------------------------
28 #ifndef THTM_FluidsH
29 #define THTM_FluidsH
30 // ---------------------------------------------------------------------------
31 
32 #include "Globales.h"
33 
34 struct stHTM_Fluid {
35  double mu;
36  double rho;
37  double Cp;
38  double k;
39  double Pr;
40  double Re;
41  double AF;
42  double Hum;
43  double R;
44  double g;
45 
46  stHTM_Fluid() {
47  }
48  ;
49 
50  double Property(double T, dVector A) {
51  return A[0] + A[1] * T + A[2] * pow2(T) + A[3] * pow3(T) + A[4] * log(T) + A[5] / T + A[6] / pow2(T);
52  }
53  ;
54 
55  double fun_mu(double T);
56 
57  double fun_k(double T);
58 
59  double fun_g(double T);
60 
61  double fun_mu(double T, dVector _X);
62 
63  double fun_k(double T, dVector _X);
64 
65  double fun_g(double T, dVector _X);
66 
67  double fun_rho(double T);
68 
69  double fun_rho(double p, double T);
70 
71  double fun_Pr(double T);
72 
73  double fun_Pr(double T, dVector _X);
74 
75  void virtual CalcProperties(double p, double T) = 0;
76 };
77 
79 
80  stHTMair() :
81  stHTM_Fluid() {
82  AF = 0.;
83  Hum = 0.;
84  }
85 
86  double fun_mu(double T) {
87  dVector A;
88  A.resize(7);
89  A[0] = 0.0000997217;
90  A[1] = 0.000000214683;
91  A[2] = -0.000000000226615;
92  A[3] = 1.1354e-13;
93  A[4] = -0.0000218018;
94  A[5] = -0.001168773;
95  A[6] = 0;
96 
97  return Property(T, A);
98  }
99 
100  double fun_k(double T) {
101  dVector A;
102  A.resize(7);
103  A[0] = -0.028882591;
104  A[1] = 0.0000746349;
105  A[2] = -0.000000020381;
106  A[3] = 6.66102E-13;
107  A[4] = 0.005818895;
108  A[5] = 0.407241757;
109  A[6] = 0;
110 
111  return Property(T, A);
112  }
113 
114  double fun_Cp(double T) {
115 
116  double Cp, x;
117 
118  if(AF == 0.0) {
119  x = Hum / (1 + Hum);
120  Cp = (1 - x) * fun_Cp_AirS(T) + x * fun_Cp_W(T);
121  } else {
122  x = 1 / AF;
123  Cp = (1 - x) * fun_Cp_AirS(T) + x * fun_Cp_Gas(T);
124  }
125  return Cp;
126  }
127 
128  double fun_Cp_AirS(double T) {
129  return -10.4199 * sqrt(T) + 2809.87 - 67227.1 / sqrt(T) + 917124.4 / T - 4174853.6 / pow150(T);
130  }
131 
132  double fun_Cp_W(double T) {
133  return -41.9055 * sqrt(T) + 10447.493 - 382002.49 / sqrt(T) + 6456647.7 / T - 37951136.5 / pow150(T);
134  }
135 
136  double fun_Cp_Gas(double T) {
137  return 926.554 + 0.43045 * T - 0.0001125 * pow2(T) + 0.000000008979 * pow3(T);
138  }
139 
140  double fun_g(double T) {
141  return fun_Cp(T) / (fun_Cp(T) - fun_R());
142  }
143 
144  double fun_R() {
145  double x, R_Air;
146 
147  if(AF == 0.0) {
148  x = Hum / (1 + Hum);
149  R_Air = (1 - x) * 287 + x * 462.1762;
150  } else {
151  x = 1 / AF;
152  R_Air = (1 - x) * 287 + x * 285.4;
153  }
154  return R_Air;
155  }
156 
157  double fun_rho(double p, double T) {
158  return (100000 * p) / (fun_R() * T);
159  }
160 
161  double fun_Pr(double T) {
162  return fun_mu(T) * fun_Cp(T) / fun_k(T);
163  }
164 
165  void CalcProperties(double p, double T) {
166  mu = fun_mu(T);
167  R = fun_R();
168  Cp = fun_Cp(T);
169  g = fun_g(T);
170  rho = fun_rho(p, T);
171  k = fun_k(T);
172  Pr = mu * Cp / k;
173  }
174 };
175 
177 
178  double mu_c1;
179  double mu_c2;
180  double mu_c3;
181 
182  stHTMoil() :
183  stHTM_Fluid() {
184  mu_c1 = 0.0115374825;
185  mu_c2 = 62.1420218;
186  mu_c3 = 275.888115;
187  }
188 
189  double fun_rho(double T) {
190  dVector A;
191  A.resize(7);
192  A[0] = -16572.74497;
193  A[1] = -38.71871842;
194  A[2] = 0.054641551;
195  A[3] = -3.45327E-05;
196  A[4] = 4398.316892;
197  A[5] = 0;
198  A[6] = 0;
199 
200  return Property(T, A);
201  }
202 
203  double fun_Cp(double T) {
204  dVector A;
205  A.resize(7);
206  A[0] = -34040.12013;
207  A[1] = -71.24978467;
208  A[2] = 0.107361063;
209  A[3] = -6.63966E-05;
210  A[4] = 8670.542844;
211  A[5] = 0;
212  A[6] = 0;
213 
214  return Property(T, A);
215  }
216 
217  double fun_mu(double T) {
218  return mu_c1 * exp(mu_c2 / (T - mu_c3));
219  }
220 
221  double fun_k(double T) {
222  dVector A;
223  A.resize(7);
224  A[0] = -13.78897843;
225  A[1] = -0.028315003;
226  A[2] = 3.80794E-05;
227  A[3] = -2.25321E-08;
228  A[4] = 3.437962885;
229  A[5] = 0;
230  A[6] = 0;
231 
232  return Property(T, A);
233  }
234 
235  double fun_Pr(double T) {
236  return fun_mu(T) * fun_Cp(T) / fun_k(T);
237  }
238 
239  void CalcProperties(double p, double T) {
240  mu = fun_mu(T);
241  Cp = fun_Cp(T);
242  rho = fun_rho(T);
243  k = fun_k(T);
244  Pr = mu * Cp / k;
245  }
246 
247 };
248 
250 
251  stHTMwater() :
252  stHTM_Fluid() {
253  }
254 
255  double fun_mu(double T) {
256 
257  return 3.117996E-11 * pow2(pow2(T)) - 4.275045E-08 * pow3(T) + 2.202288E-05 * pow2(T) - 5.057924E-03 * T + 4.378690E-01;
258  }
259 
260  double fun_Cp(double T) {
261 // return -41.9055 * sqrt(T) + 10447.493 - 382002.49 / sqrt(T)
262 // + 6456647.7 / T - 37951136.5 / pow150(T);
263  double T_K = __units::degCToK(T);
264  return 1987087.20115782 + T_K * (4493.2068586048 + T_K * (-6.72837895716355 + T_K * 0.00446936120713534)) -
265  498966.339389734 * log(T_K);
266  }
267 
268  double fun_rho(double p, double T) {
269  return -1.374183E-07 * pow2(pow2(T)) + 1.937996E-04 * pow3(T) - 1.050544E-01 * pow2(
270  T) + 2.527211E+01 * T - 1.249635E+03;
271  }
272 
273  double fun_k(double T) {
274  return 0.58;
275  }
276 
277  double fun_Pr(double T) {
278  return fun_mu(T) * fun_Cp(T) / fun_k(T);
279  }
280 
281  void CalcProperties(double p, double T) {
282  mu = fun_mu(T);
283  Cp = fun_Cp(T);
284  rho = fun_rho(p, T);
285  k = fun_k(T);
286  Pr = mu * Cp / k;
287  }
288 
289 };
290 
292 
293  stHTMFreshair() :
294  stHTM_Fluid() {
295  }
296 
297  double fun_mu(double T) {
298  dVector A;
299  A.resize(7);
300  A[0] = 0.0000997217;
301  A[1] = 0.000000214683;
302  A[2] = -0.000000000226615;
303  A[3] = 1.1354e-13;
304  A[4] = -0.0000218018;
305  A[5] = -0.001168773;
306  A[6] = 0;
307 
308  return Property(T, A);
309  }
310 
311  double fun_k(double T) {
312  dVector A;
313  A.resize(7);
314  A[0] = -0.028882591;
315  A[1] = 0.0000746349;
316  A[2] = -0.000000020381;
317  A[3] = 6.66102E-13;
318  A[4] = 0.005818895;
319  A[5] = 0.407241757;
320  A[6] = 0;
321 
322  return Property(T, A);
323  }
324 
325  double fun_Cp(double T) {
326  return -10.4199 * sqrt(T) + 2809.87 - 67227.1 / sqrt(T) + 917124.4 / T - 4174853.6 / pow150(T);
327  }
328 
329  double fun_g(double T) {
330  return fun_Cp(T) / (fun_Cp(T) - fun_R());
331  }
332 
333  double fun_R() {
334  double x, R_Air;
335 
336  R_Air = 287.;
337 
338  return R_Air;
339  }
340 
341  double fun_rho(double p, double T) {
342  return (100000 * p) / (fun_R() * T);
343  }
344 
345  double fun_Pr(double T) {
346  return fun_mu(T) * fun_Cp(T) / fun_k(T);
347  }
348 
349  void CalcProperties(double p, double T) {
350  mu = fun_mu(T);
351  R = fun_R();
352  Cp = fun_Cp(T);
353  g = fun_g(T);
354  rho = fun_rho(p, T);
355  k = fun_k(T);
356  Pr = mu * Cp / k;
357  }
358 };
359 
361 
362  stHTMBurntGas() :
363  stHTM_Fluid() {
364  }
365 
366  double fun_mu(double T) {
367  dVector A;
368  A.resize(7);
369  A[0] = 0.0000997217;
370  A[1] = 0.000000214683;
371  A[2] = -0.000000000226615;
372  A[3] = 1.1354e-13;
373  A[4] = -0.0000218018;
374  A[5] = -0.001168773;
375  A[6] = 0;
376 
377  return Property(T, A);
378  }
379 
380  double fun_k(double T) {
381  dVector A;
382  A.resize(7);
383  A[0] = -0.028882591;
384  A[1] = 0.0000746349;
385  A[2] = -0.000000020381;
386  A[3] = 6.66102E-13;
387  A[4] = 0.005818895;
388  A[5] = 0.407241757;
389  A[6] = 0;
390 
391  return Property(T, A);
392  }
393 
394  double fun_Cp(double T) {
395  return 926.554 + 0.43045 * T - 0.0001125 * pow2(T) + 0.000000008979 * pow3(T);
396  }
397 
398  double fun_g(double T) {
399  return fun_Cp(T) / (fun_Cp(T) - fun_R());
400  }
401 
402  double fun_R() {
403  double x, R_Air;
404 
405  R_Air = 285.4;
406 
407  return R_Air;
408  }
409 
410  double fun_rho(double p, double T) {
411  return (100000 * p) / (fun_R() * T);
412  }
413 
414  double fun_Pr(double T) {
415  return fun_mu(T) * fun_Cp(T) / fun_k(T);
416  }
417 
418  void CalcProperties(double p, double T) {
419  mu = fun_mu(T);
420  R = fun_R();
421  Cp = fun_Cp(T);
422  g = fun_g(T);
423  rho = fun_rho(p, T);
424  k = fun_k(T);
425  Pr = mu * Cp / k;
426  }
427 };
428 
430 
431  dVector X;
432  stHTMwater X0_wat;
433  stHTMBurntGas X1_gas;
434  stHTMFreshair X2_air;
435 
436  stBurntAirWater() :
437  stHTM_Fluid() {
438  X.resize(2);
439  X[0] = 0.;
440  X[1] = 0.;
441  }
442 
444  stHTM_Fluid() {
445  X.resize(2);
446  X[0] = _X[0];
447  X[1] = _X[1];
448  }
449 
450  void Put_comp(dVector _X) {
451 
452  X[0] = _X[0];
453  X[1] = _X[1];
454  }
455 
456  double fun_mu(double T) {
457  return X[0] * X0_wat.fun_mu(T) + X[1] * X1_gas.fun_mu(T) + (1 - X[0] - X[1]) * X2_air.fun_mu(T);
458  }
459 
460  double fun_mu(double T, dVector _X) {
461  return _X[0] * X0_wat.fun_mu(T) + _X[1] * X1_gas.fun_mu(T) + (1 - X[0] - X[1]) * X2_air.fun_mu(T);
462  }
463 
464  double fun_Cp(double T) {
465  return X[0] * X0_wat.fun_Cp(T) + X[1] * X1_gas.fun_Cp(T) + (1 - X[0] - X[1]) * X2_air.fun_Cp(T);
466  }
467 
468  double fun_Cp(double T, dVector _X) {
469  return _X[0] * X0_wat.fun_Cp(T) + _X[1] * X1_gas.fun_Cp(T) + (1 - X[0] - X[1]) * X2_air.fun_Cp(T);
470  }
471 
472  double fun_R() {
473  return X[1] * X1_gas.fun_R() + (1 - X[1]) * X2_air.fun_R();
474  }
475 
476  double fun_R(dVector _X) {
477  return _X[1] * X1_gas.fun_R() + (1 - _X[1]) * X2_air.fun_R();
478  }
479 
480  double fun_Pr(double T) {
481  return fun_mu(T) * fun_Cp(T) / fun_k(T);
482  }
483 
484  double fun_Pr(double T, dVector _X) {
485  return fun_mu(T, _X) * fun_Cp(T, _X) / fun_k(T, _X);
486  }
487 
488  double fun_g(double T) {
489  return fun_Cp(T) / (fun_Cp(T) - fun_R());
490  }
491 
492  double fun_g(double T, dVector _X) {
493  return fun_Cp(T, _X) / (fun_Cp(T, _X) - fun_R(_X));
494  }
495 
496  void CalcProperties(double p, double T) {
497 
498  }
499 
500 };
501 
502 #endif
503 
pow150
T pow150(T x)
Returns x to the power of 1.5.
Definition: Math_wam.h:140
stHTMair
Definition: THTM_Fluids.h:78
stHTMwater
Definition: THTM_Fluids.h:249
pow3
T pow3(T x)
Returns x to the power of 3.
Definition: Math_wam.h:101
stHTMoil
Definition: THTM_Fluids.h:176
stHTMFreshair
Definition: THTM_Fluids.h:291
stHTM_Fluid
Definition: THTM_Fluids.h:34
pow2
T pow2(T x)
Returns x to the power of 2.
Definition: Math_wam.h:88
stHTMBurntGas
Definition: THTM_Fluids.h:360
dVector
std::vector< double > dVector
Double vector.
Definition: Math_wam.h:70
stBurntAirWater
Definition: THTM_Fluids.h:429