1 //////////////////////////////////////////////////////////////////////////////
\r
7 //////////////////////////////////////////////////////////////////////////////
\r
9 #ifndef fluid_included
\r
10 #define fluid_included
\r
14 #include <SFichier.h>
\r
15 #include <Process.h>
\r
20 fluid() { fluide = 0; };
\r
21 void set_fluid(Nom &nom)
\r
23 Noms liste(2); liste(0) = "Na"; liste(1) = "LBe";
\r
24 fluide = liste.rang(nom);
\r
27 if (!Process::me()) Cerr << "fluide " << nom << " non reconnu parmi : " << liste << finl;
\r
32 //proprietes physiques
\r
33 /* Lois physiques du sodium issues de fits sur DEN/DANS/DM2S/STMF/LMES/RT/12-018/A */
\r
35 #define Tk ((T) + 273.15)
\r
36 #define Ksil (fluide ? 3.022e-11 : 1e-9)
\r
38 /* inverse de la densite du liquide (hors pression) */
\r
39 inline double IRhoL(double T)
\r
41 return fluide ? (9.03e-5 + Tk * (1.003e-8 + Tk * 2.01e-12))
\r
42 : (9.829728e-4 + Tk * (2.641186e-7 + Tk * (-3.340743e-11 + Tk * 6.680973e-14)));
\r
46 inline double DTIRhoL(double T)
\r
48 return fluide ? (1.003e-8 + 2.01e-12 * 2)
\r
49 : (2.641186e-7 + Tk * (-3.340743e-11 * 2 + Tk * 6.680973e-14 * 3));
\r
52 /* densite du liquide */
\r
53 inline double RhoL(double T, double P)
\r
55 return exp(Ksil * (P - 1e5)) / IRhoL(T);
\r
58 /* enthalpie du liquide */
\r
59 inline double HL(double T, double P)
\r
61 return (fluide ? 97.98e3 + Tk * (159 + Tk * (-2.72e-2 / 2 + Tk * 7.12e-6 / 3))
\r
62 : 2992600 / Tk - 365487.2 + Tk * (1658.2 + Tk * (-.42395 + Tk * 1.4847e-4)))
\r
63 + (IRhoL(T) - Tk * DTIRhoL(T)) * (1 - exp(Ksil*(1e5 - P))) / Ksil;
\r
67 inline double DTHL(double T, double P)
\r
69 return (fluide ? 159 + Tk * (-2.72e-2 + Tk * 7.12e-6)
\r
70 : -2992600 / (Tk*Tk) + 1658.2 + Tk * (-.42395 * 2 + Tk * 1.4847e-4 * 3))
\r
71 - Tk * (fluide ? 0 : -3.340743e-11 * 2 + Tk * 6.680973e-14 * 6) * (1 - exp(Ksil*(1e5 - P))) / Ksil;
\r
74 /* inverses par methode de Newton */
\r
75 inline double TLh(double h, double T0, double P)
\r
77 double T = T0, sec; int i;
\r
78 for (i = 0; i < 100 && dabs( sec = HL(T, P) - h ) > 1e-8; i++)
\r
79 T -= sec / DTHL(T, P);
\r
80 return i < 100 ? T : -1e10;
\r
83 inline double DTRhoL(double T, double P)
\r
85 return - exp(Ksil * (P - 1e5)) * DTIRhoL(T) / (IRhoL(T) * IRhoL(T));
\r
87 inline double DPRhoL(double T, double P)
\r
89 return Ksil * RhoL(T, P);
\r
93 inline double DPHL(double T, double P)
\r
95 return (IRhoL(T) - Tk * DTIRhoL(T)) * exp(Ksil * (1e5 - P));
\r
98 /* energie volumique du liquide */
\r
99 inline double EL(double T, double P)
\r
101 return HL(T, P) * RhoL(T, P);
\r
104 inline double DTEL(double T, double P)
\r
106 return DTHL(T, P) * RhoL(T, P) + HL(T, P) * DTRhoL(T, P);
\r
108 inline double DPEL(double T, double P)
\r
110 return DPHL(T, P) * RhoL(T, P) + HL(T, P) * DPRhoL(T, P);
\r
113 inline double TLe(double e, double T0, double P)
\r
115 double T = T0, sec; int i;
\r
116 for (i = 0; i < 100 && dabs( sec = EL(T, P) - e ) > 1e-3; i++)
\r
117 T -= sec / DTEL(T, P);
\r
118 return i < 100 ? T : -1e10;
\r
121 /* viscosite du liquide*/
\r
122 inline double MuL( double T )
\r
124 return fluide ? 4.94e-4 * exp(754.1 / Tk)
\r
125 : exp( -6.4406 - 0.3958 * log(Tk) + 556.835 / Tk );
\r
127 /* conductivite du liquide */
\r
128 inline double LambdaL( double T )
\r
130 return fluide ? 3.61 + Tk * (1.517e-2 - Tk * 1.741e-6)
\r
131 : max(124.67 + Tk * (-.11381 + Tk * (5.5226e-5 - Tk * 1.1848e-8)), 0.045);
\r
133 /* tension superficielle */
\r
134 inline double SigmaL(double T)
\r
136 return 0.2405 * pow(1 - Tk / 2503.7, 1.126);
\r
138 /* Nusselt du liquide -> Skupinski */
\r
139 inline double NuL(double T, double P, double w, double Dh)
\r
141 double Pe = RhoL(T, P) * dabs(w) * Dh * DTHL(T, P) / LambdaL(T);
\r
142 return 4.82 + 0.0185 * pow(Pe, 0.827);
\r
144 /* temperature de saturation */
\r
145 inline double Tsat( double P )
\r
147 double A = 7.8130, B = 11209, C = 5.249e5;
\r
148 return 2 * C / ( -B + sqrt(B*B + 4 * A * C - 4 * C * log(P / 1e6))) - 273.15;
\r
151 inline double DTsat( double P )
\r
153 double A = 7.8130, B = 11209, C = 5.249e5, Tsk = Tsat(P) + 273.15;
\r
154 return Tsk * Tsk / ( P * sqrt(B*B + 4 * A * C - 4 * C * log(P / 1e6)));
\r
156 /* pression de vapeur saturante */
\r
157 inline double Psat( double T )
\r
159 double A = 7.8130, B = 11209, C = 5.249e5;
\r
160 return 1e6 * exp(A - B / Tk - C / (Tk * Tk));
\r
163 inline double DPsat( double T )
\r
165 double B = 11209, C = 5.249e5;
\r
166 return (B / (Tk * Tk) + 2 * C / (Tk * Tk * Tk)) * Psat(T);
\r
168 /* enthalpie massique de saturation */
\r
169 inline double Hsat( double P )
\r
171 return HL(Tsat(P), P);
\r
174 inline double DHsat( double P )
\r
176 return DTsat(P) * DTHL(Tsat(P), P) + DPHL(Tsat(P), P);
\r
178 /* chaleur latente, a prendre a Tsat */
\r
179 inline double Lvap( double P )
\r
181 double Tc = 2503.7, Tsk = Tsat(P) + 273.15;
\r
182 return 3.9337e5 * ( 1 - Tsk / Tc) + 4.3986e6 * pow( 1 - Tsk / Tc, .29302);
\r
186 inline double DLvap( double P )
\r
188 double Tc = 2503.7, Tsk = Tsat(P) + 273.15;
\r
189 return DTsat(P) * (-3.9337e5 / Tc - 4.3986e6 * .29302 * pow( 1 - Tsk / Tc, .29302 - 1) / Tc);
\r
192 /* densite de la vapeur : (gaz parfait) * f1(Tsat) * f2(DTsat)*/
\r
193 #define f1(Ts) (2.49121 + Ts * (-5.53796e-3 + Ts * (7.5465e-6 + Ts * (-4.20217e-9 + Ts * 8.59212e-13))))
\r
194 #define Df1(Ts) (-5.53796e-3 + Ts * (2 * 7.5465e-6 + Ts * (-3 * 4.20217e-9 + Ts * 4 * 8.59212e-13)))
\r
195 #define f2(dT) (1 + dT * (-5e-4 + dT * (6.25e-7 - dT * 4.1359e-25)))
\r
196 #define Df2(dT) (-5e-4 + dT * (2 * 6.25e-7 - dT * 3 * 4.1359e-25))
\r
197 inline double RhoV( double T, double P )
\r
199 double Tsk = Tsat(P) + 273.15, dT = Tk - Tsk;
\r
200 return P * 2.7650313e-3 * f1(Tsk) * f2(dT) / Tk;
\r
203 inline double DTRhoV( double T, double P )
\r
205 double Tsk = Tsat(P) + 273.15, dT = Tk - Tsk;
\r
206 return P * 2.7650313e-3 * f1(Tsk) * (Df2(dT) - f2(dT) / Tk) / Tk;
\r
208 inline double DPRhoV( double T, double P )
\r
210 double Tsk = Tsat(P) + 273.15, dT = Tk - Tsk;
\r
211 return 2.7650313e-3 * (f1(Tsk) * f2(dT) + P * DTsat(P) * (Df1(Tsk) * f2(dT) - f1(Tsk) * Df2(dT))) / Tk;
\r
217 /* inverse de la densite de la vapeur */
\r
218 inline double IRhoV( double T, double P )
\r
220 return 1 / RhoV(T, P);
\r
223 inline double DTIRhoV( double T, double P )
\r
225 return -DTRhoV(T, P) / (RhoV(T, P) * RhoV(T, P));
\r
227 inline double DPIRhoV( double T, double P )
\r
229 return -DPRhoV(T, P) / (RhoV(T, P) * RhoV(T, P));
\r
231 /* enthalpie de la vapeur */
\r
232 #define CpVs(Ts) (-1223.89 + Ts * (14.1191 + Ts * (-1.62025e-2 + Ts * 5.76923e-6)))
\r
233 #define DCpVs(Ts) (14.1191 + Ts * (- 2 * 1.62025e-2 + Ts * 3 * 5.76923e-6))
\r
234 inline double HV( double T, double P)
\r
236 double Ts = Tsat(P), dT = T - Ts, Cp0 = 910, k = -4.6e-3;
\r
237 return HL(Ts, P) + Lvap(P) + Cp0 * dT + (Cp0 - CpVs(Ts)) * (1 - exp(k * dT)) / k;
\r
240 inline double DTHV( double T, double P)
\r
242 double Ts = Tsat(P), dT = T - Ts, Cp0 = 910, k = -4.6e-3;
\r
243 return Cp0 - (Cp0 - CpVs(Ts)) * exp(k * dT);
\r
245 inline double DPHV( double T, double P)
\r
247 double Ts = Tsat(P), dT = T - Ts, Cp0 = 910, k = -4.6e-3;
\r
248 return DPHL(Ts, P) + DLvap(P) + DTsat(P) * (DTHL(Ts, P) - Cp0 - DCpVs(Ts) * (1 - exp(k * dT)) / k + (Cp0 - CpVs(Ts)) * exp(k * dT));
\r
252 /* energie volumique de la vapeur */
\r
253 inline double EV( double T, double P)
\r
255 return RhoV(T, P) * HV(T, P);
\r
258 inline double DTEV(double T, double P)
\r
260 return DTRhoV(T, P) * HV(T, P) + RhoV(T, P) * DTHV(T, P);
\r
262 inline double DPEV(double T, double P)
\r
264 return DPRhoV(T, P) * HV(T, P) + RhoV(T, P) * DPHV(T, P);
\r
267 /* conductivite de la vapeur */
\r
268 inline double LambdaV( double T )
\r
273 /* viscosite de la vapeur */
\r
274 inline double MuV( double T )
\r
276 return 2.5e-6 + 1.5e-8 * Tk;
\r
279 /* Nusselt de la vapeur -> Dittus/ Boetler */
\r
280 inline double NuV(double T, double P, double w, double Dh)
\r
282 double Re = RhoV(T, P) * dabs(w) * Dh / MuV(T), Pr = MuV(T) * DTHV(T, P) / LambdaV(T);
\r
283 return 0.023 * pow(Re, 0.8) * pow(Pr, 0.4);
\r
286 inline double Rho( double T, double x, double P )
\r
288 return 1 / ( (1-x)/RhoL(T, P) + x / RhoV(T, P) );
\r