TrioCFD 1.9.8
TrioCFD documentation
Loading...
Searching...
No Matches
Lois_sodium.h
1/****************************************************************************
2* Copyright (c) 2023, CEA
3* All rights reserved.
4*
5* Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
6* 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
7* 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
8* 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
9*
10* THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
11* IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
12* OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
13*
14*****************************************************************************/
15
16#ifndef lois_sodium_span_included
17#define lois_sodium_span_included
18
19#include <assert.h>
20#include <span.hpp>
21#include <vector>
22#include <cmath>
23
24/* Lois physiques du sodium issues de fits sur DEN/DANS/DM2S/STMF/LMES/RT/12-018/A */
25
26// iterator index !
27#define i_it std::distance(res.begin(), &val)
28
29using VectorD = std::vector<double>;
30using SpanD = tcb::span<double>;
31
32const double Ksil = 1e-9;
33
34/* inverse de la densite du liquide (hors pression) */
35inline void IRhoL(const SpanD T, SpanD res, int ncomp, int id)
36{
37 for (auto& val : res)
38 {
39 const double tk = T[i_it * ncomp + id] + 273.15;
40 val = (9.829728e-4 + tk * (2.641186e-7 + tk * (-3.340743e-11 + tk * 6.680973e-14)));
41 }
42}
43
44/* sa derivee */
45inline void DTIRhoL(const SpanD T, SpanD res, int ncomp, int id)
46{
47 for (auto& val : res)
48 {
49 const double tk = T[i_it * ncomp + id] + 273.15;
50 val = (2.641186e-7 + tk * (-3.340743e-11 * 2 + tk * 6.680973e-14 * 3));
51 }
52}
53
54inline void RhoL(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
55{
56 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
57 IRhoL(T,res,ncomp,id); // Attention on utilise res ...
58 for (auto& val : res) val = exp(Ksil * (P[i_it] - 1e5)) / val;
59}
60
61/* derivees */
62inline void DPRhoL(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
63{
64 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
65 RhoL(T,P,res,ncomp,id); // Attention on utilise res ...
66 for (auto& val : res) val = Ksil * val;
67}
68
69inline void DTRhoL(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
70{
71 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
72 VectorD invRhoL((int)res.size()), DinvRhoL((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
73 IRhoL(T,SpanD(invRhoL),ncomp,id);
74 DTIRhoL(T,SpanD(DinvRhoL),ncomp,id);
75 for (auto& val : res) val = - exp(Ksil * (P[i_it] - 1e5)) * DinvRhoL[i_it] / (invRhoL[i_it] * invRhoL[i_it]);
76}
77
78/* enthalpie du liquide */
79inline void HL(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
80{
81 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
82 VectorD invRhoL((int)res.size()), DinvRhoL((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
83 IRhoL(T,SpanD(invRhoL),ncomp,id);
84 DTIRhoL(T,SpanD(DinvRhoL),ncomp,id);
85 for (auto& val : res)
86 {
87 const double tk = T[i_it * ncomp + id] + 273.15;
88 val = (2992600 / tk - 365487.2 + tk * (1658.2 + tk * (-.42395 + tk * 1.4847e-4)))+ (invRhoL[i_it] - tk * DinvRhoL[i_it]) * (1 - exp(Ksil*(1e5 - P[i_it]))) / Ksil;
89 }
90}
91
92/* derivees */
93inline void DPHL(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
94{
95 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
96 VectorD invRhoL((int)res.size()), DinvRhoL((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
97 IRhoL(T,SpanD(invRhoL),ncomp,id);
98 DTIRhoL(T,SpanD(DinvRhoL),ncomp,id);
99 for (auto& val : res)
100 {
101 const double tk = T[i_it * ncomp + id] + 273.15;
102 val = (invRhoL[i_it] - tk * DinvRhoL[i_it]) * exp(Ksil * (1e5 - P[i_it]));
103 }
104}
105
106inline void DTHL(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
107{
108 assert ((int)T.size() == ncomp * (int)res.size() && (int)T.size() == ncomp * (int)P.size());
109 for (auto& val : res)
110 {
111 const double tk = T[i_it * ncomp + id] + 273.15;
112 val = (-2992600 / (tk*tk) + 1658.2 + tk * (-.42395 * 2 + tk * 1.4847e-4 * 3))- tk * (-3.340743e-11 * 2 + tk * 6.680973e-14 * 6) * (1 - exp(Ksil*(1e5 - P[i_it]))) / Ksil;
113 }
114}
115
116inline void BETAL(const SpanD T, SpanD res, int ncomp, int id)
117{
118 assert ((int)T.size() == ncomp * (int)res.size());
119 VectorD invRhoL((int)res.size()), DinvRhoL((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
120 IRhoL(T,SpanD(invRhoL),ncomp,id);
121 DTIRhoL(T,SpanD(DinvRhoL),ncomp,id);
122 for (auto& val : res) val = DinvRhoL[i_it] / invRhoL[i_it];
123}
124
125/* viscosite du liquide*/
126inline void MuL(const SpanD T, SpanD res, int ncomp, int id)
127{
128 assert ((int)T.size() == ncomp * (int)res.size());
129 for (auto& val : res)
130 {
131 const double tk = T[i_it * ncomp + id] + 273.15;
132 val = exp( -6.4406 - 0.3958 * log(tk) + 556.835 / tk );
133 }
134}
135
136/* conductivite du liquide */
137inline void LambdaL(const SpanD T, SpanD res, int ncomp, int id)
138{
139 assert ((int)T.size() == ncomp * (int)res.size());
140 for (auto& val : res)
141 {
142 const double tk = T[i_it * ncomp + id] + 273.15;
143 val = std::max(124.67 + tk * (-.11381 + tk * (5.5226e-5 - tk * 1.1848e-8)), 0.045);
144 }
145}
146
147/* tension superficielle */
148inline void SigmaL(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
149{
150 assert ((int)P.size() == (int)res.size());
151 for (auto& val : res)
152 {
153 const double tk = T[i_it * ncomp + id] + 273.15;
154 val = 0.2405 * pow(1 - tk / 2503.7, 1.126);
155 }
156}
157
158/*
159 * ****************************
160 * ********** GAZ ***********
161 * ****************************
162 */
163
164/* temperature de saturation */
165inline void Tsat_Na(const SpanD P, SpanD res, int ncomp, int ind)
166{
167 assert (ncomp * (int)P.size() == (int)res.size());
168 const double A = 7.8130, B = 11209, C = 5.249e5;
169 for (int i =0; i < (int)P.size(); i++) res[i * ncomp + ind] = 2 * C / ( -B + sqrt(B*B + 4 * A * C - 4 * C * log(P[i] / 1e6))) - 273.15;
170}
171
172/* sa derivee */
173inline void DTsat_Na(const SpanD P, SpanD res, int ncomp, int ind)
174{
175 assert (ncomp * (int)P.size() == (int)res.size());
176 const double A = 7.8130, B = 11209, C = 5.249e5;
177 Tsat_Na(P,res,ncomp,ind); // XXX : si c'est complique va falloir creer un tableau local ...
178 for (int i =0; i < (int)P.size(); i++)
179 {
180 const double Tsk = res[i * ncomp + ind] + 273.15; // XXX : need in Kelvin
181 res[i * ncomp + ind] = Tsk * Tsk / ( P[i] * sqrt(B*B + 4 * A * C - 4 * C * log(P[i] / 1e6)));
182 }
183}
184
185/* pression de vapeur saturante */
186inline void Psat_Na(const SpanD T, SpanD res, int ncomp, int ind)
187{
188 assert ((int)T.size() == (int)res.size());
189 const double A = 7.8130, B = 11209, C = 5.249e5;
190 for (int i =0; i < (int)T.size() / ncomp; i++)
191 {
192 const double tk = T[i * ncomp + ind] + 273.15;
193 res[i * ncomp + ind] = 1.e6 * exp(A - B / tk - C / (tk * tk));
194 }
195}
196
197/* sa derivee */
198inline void DPsat_Na(const SpanD T, SpanD res, int ncomp, int ind)
199{
200 assert ((int)T.size() == (int)res.size());
201 const double B = 11209, C = 5.249e5;
202 // attention on utilise resu !
203 Psat_Na(T,res,ncomp,ind);
204 for (int i =0; i < (int)T.size() / ncomp; i++)
205 {
206 const double tk = T[i * ncomp + ind] + 273.15;
207 res[i * ncomp + ind] = (B / (tk * tk) + 2 * C / (tk * tk * tk)) * res[i * ncomp + ind];
208 }
209}
210
211/* enthalpie massique de saturation */
212inline void Hsat(const SpanD P, SpanD res, int ncomp, int ind)
213{
214 assert (ncomp * (int)P.size() == (int)res.size());
215 VectorD Ts_((int)P.size()), HL_((int)P.size());
216 Tsat_Na(P,SpanD(Ts_),1,0); // 1D
217 HL(SpanD(Ts_),P,SpanD(HL_),1,0); // 1D
218 for (int i =0; i < (int)P.size(); i++) res[i * ncomp + ind] = HL_[i];
219}
220
221/* sa derivee */
222inline void DHsat(const SpanD P, SpanD res, int ncomp, int ind)
223{
224 assert (ncomp * (int)P.size() == (int)res.size());
225 VectorD Ts_((int)P.size()), DTs_((int)P.size()), DTHL_((int)P.size()), DPHL_((int)P.size());
226 Tsat_Na(P,SpanD(Ts_),1,0); // 1D
227 DTsat_Na(P,SpanD(DTs_),1,0); // 1D
228 DTHL(SpanD(Ts_),P,SpanD(DTHL_),1,0);
229 DPHL(SpanD(Ts_),P,SpanD(DPHL_),1,0);
230 for (int i =0; i < (int)P.size(); i++) res[i * ncomp + ind] = DTs_[i] * DTHL_[i] + DPHL_[i];
231}
232
233/* densite de la vapeur : (gaz parfait) * f1(Tsat_Na) * f2(DTsat_Na)*/
234#define f1 (2.49121 + Tsk * (-5.53796e-3 + Tsk * (7.5465e-6 + Tsk * (-4.20217e-9 + Tsk * 8.59212e-13))))
235#define f2 (1 + dT * (-5e-4 + dT * (6.25e-7 - dT * 4.1359e-25)))
236#define Df1 (-5.53796e-3 + Tsk * (2 * 7.5465e-6 + Tsk * (-3 * 4.20217e-9 + Tsk * 4 * 8.59212e-13)))
237#define Df2 (-5e-4 + dT * (2 * 6.25e-7 - dT * 3 * 4.1359e-25))
238
239inline void RhoV(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
240{
241 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
242 Tsat_Na(P,res,1,0); // XXX : attention on utilise res
243 for (auto& val : res)
244 {
245 const double Tsk = val + 273.15; // XXX : need in Kelvin
246 const double tk = T[i_it * ncomp + id] + 273.15;
247 const double dT = tk - Tsk;
248 val = P[i_it] * 2.7650313e-3 * f1 * f2 / tk;
249 }
250}
251
252/* ses derivees */
253inline void DPRhoV(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
254{
255 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
256 VectorD Tsk_((int)res.size()), DTsk_((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
257 Tsat_Na(P,SpanD(Tsk_),1,0);
258 DTsat_Na(P,SpanD(DTsk_),1,0);
259 for (auto& val : res)
260 {
261 const double Tsk = Tsk_[i_it] + 273.15; // XXX : need in Kelvin
262 const double tk = T[i_it * ncomp + id] + 273.15;
263 const double dT = tk - Tsk;
264 val = 2.7650313e-3 * (f1 * f2 + P[i_it] * DTsk_[i_it] * (Df1 * f2 - f1 * Df2)) / tk;
265 }
266}
267
268inline void DTRhoV(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
269{
270 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
271 Tsat_Na(P,res,1,0); // XXX : attention on utilise res
272 for (auto& val : res)
273 {
274 const double Tsk = val + 273.15; // XXX : need in Kelvin
275 const double tk = T[i_it * ncomp + id] + 273.15;
276 const double dT = tk - Tsk;
277 val = P[i_it] * 2.7650313e-3 * f1 * (Df2 - f2 / tk) / tk;
278 }
279}
280
281#undef f1
282#undef Df1
283#undef f2
284#undef Df2
285
286/* chaleur latente, a prendre a Tsat_Na */
287inline void Lvap_Na(const SpanD P, SpanD res, int ncomp, int ind)
288{
289 assert (ncomp * (int)P.size() == (int)res.size());
290 Tsat_Na(P,res,ncomp,ind); // XXX : attention on utilise res. XXX : si c'est complique va falloir creer un tableau local ...
291 const double Tc = 2503.7;
292 for (int i =0; i < (int)P.size(); i++)
293 {
294 const double Tsk = res[i * ncomp + ind] + 273.15; // XXX : need in Kelvin
295 res[i * ncomp + ind] = 3.9337e5 * ( 1 - Tsk / Tc) + 4.3986e6 * pow( 1 - Tsk / Tc, .29302);
296 }
297}
298
299/* sa derivee */
300inline void DLvap_Na(const SpanD P, SpanD res, int ncomp, int ind)
301{
302 assert (ncomp * (int)P.size() == (int)res.size());
303 const double Tc = 2503.7;
304 VectorD Tsk_((int)P.size()), DTsk_((int)P.size());
305 Tsat_Na(P,SpanD(Tsk_),1,0); // 1D
306 DTsat_Na(P,SpanD(DTsk_),1,0); // 1D
307 for (int i =0; i < (int)P.size(); i++)
308 {
309 const double Tsk = Tsk_[i] + 273.15; // XXX : need in Kelvin
310 res[i * ncomp + ind] = DTsk_[i] * (-3.9337e5 / Tc - 4.3986e6 * .29302 * pow( 1 - Tsk / Tc, .29302 - 1) / Tc);
311 }
312}
313
314/* enthalpie massique de saturation (vapeur = Hsat(P) + Lvap_Na(P)) */
315inline void HVsat(const SpanD P, SpanD res, int ncomp, int ind)
316{
317 assert (ncomp * (int)P.size() == (int)res.size());
318 VectorD Hs_((int)P.size()), Lv_((int)P.size());
319 Hsat(P,SpanD(Hs_),1,0);
320 Lvap_Na(P,SpanD(Lv_),1,0);
321 for (int i =0; i < (int)P.size(); i++) res[i * ncomp + ind] = Hs_[i] + Lv_[i];
322}
323
324/* sa derivee (vapeur = DHsat(P) + DLvap_Na(P)) */
325inline void DHVsat(const SpanD P, SpanD res, int ncomp, int ind)
326{
327 assert (ncomp * (int)P.size() == (int)res.size());
328 VectorD DHs_((int)P.size()), DLv_((int)P.size());
329 DHsat(P,SpanD(DHs_),1,0);
330 DLvap_Na(P,SpanD(DLv_),1,0);
331 for (int i =0; i < (int)P.size(); i++) res[i * ncomp + ind] = DHs_[i] + DLv_[i];
332}
333
334/* enthalpie de la vapeur */
335inline void HV(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
336{
337 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
338 const double Cp0 = 910, k = -4.6e-3;
339 VectorD Ts_((int)res.size()), Lvap_((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
340 Tsat_Na(P,SpanD(Ts_),1,0);
341 Lvap_Na(P,SpanD(Lvap_),1,0);
342 HL(SpanD(Ts_),P,res,1,0); // XXX : attention on utilise res
343 for (auto& val : res)
344 {
345 const double dT = T[i_it * ncomp + id] - Ts_[i_it]; // XXX : pas meme taille attention
346 const double CpVs = (-1223.89 + Ts_[i_it] * (14.1191 + Ts_[i_it] * (-1.62025e-2 + Ts_[i_it] * 5.76923e-6)));
347 val = val + Lvap_[i_it] + Cp0 * dT + (Cp0 - CpVs) * (1 - exp(k * dT)) / k;
348 }
349}
350
351/* ses derivees */
352inline void DPHV(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
353{
354 assert ((int)T.size() == ncomp * (int)P.size() && (int)T.size() == ncomp * (int)res.size());
355 VectorD DTHL_((int)res.size()), DPHL_((int)res.size()), Ts_((int)res.size()), DTs_((int)res.size()), DLvap_((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
356 Tsat_Na(P,SpanD(Ts_),1,0);
357 DTsat_Na(P,SpanD(DTs_),1,0);
358 DLvap_Na(P,SpanD(DLvap_),1,0);
359 DTHL(SpanD(Ts_),P,SpanD(DTHL_),1,0);
360 DPHL(SpanD(Ts_),P,SpanD(DPHL_),1,0);
361
362 const double Cp0 = 910, k = -4.6e-3;
363 for (auto& val : res)
364 {
365 const double dT = T[i_it * ncomp + id] - Ts_[i_it]; // XXX : pas meme taille attention
366 const double CpVs = (-1223.89 + Ts_[i_it] * (14.1191 + Ts_[i_it] * (-1.62025e-2 + Ts_[i_it] * 5.76923e-6))); // TODO : FIXME : a faire une methode
367 const double DCpVs = (14.1191 + Ts_[i_it] * (- 2 * 1.62025e-2 + Ts_[i_it] * 3 * 5.76923e-6));
368 val = DPHL_[i_it] + DLvap_[i_it] + DTs_[i_it] * (DTHL_[i_it] - Cp0 - DCpVs * (1 - exp(k * dT)) / k + (Cp0 - CpVs) * exp(k * dT));
369 }
370}
371
372inline void DTHV(const SpanD T, const SpanD P, SpanD res, int ncomp, int id)
373{
374 assert ((int)T.size() == ncomp * (int)res.size() && (int)T.size() == ncomp * (int)P.size());
375 const double Cp0 = 910, k = -4.6e-3;
376 Tsat_Na(P,res,1,0); // XXX : attention on utilise res
377 for (auto& val : res)
378 {
379 const double dT = T[i_it * ncomp + id] - val; // XXX : pas meme taille attention : val = res[i_it]
380 const double CpVs = (-1223.89 + val * (14.1191 + val * (-1.62025e-2 + val * 5.76923e-6)));
381 val = Cp0 - (Cp0 - CpVs) * exp(k * dT);
382 }
383}
384
385/* inverse de la densite de la vapeur */
386inline void IRhoV_vec(SpanD T, const SpanD P, SpanD res, int ncomp, int id)
387{
388 RhoV(T,P,res,ncomp,id); // XXX : attention on utilise res
389 for (auto& val : res) val = 1. / val;
390}
391/* ses derivees */
392inline void DTIRhoV_vec(SpanD T, const SpanD P, SpanD res, int ncomp, int id)
393{
394 VectorD RhoV_((int)res.size()), DTRhoV_((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
395 RhoV(T,P,SpanD(RhoV_),ncomp,id);
396 DTRhoV(T,P,SpanD(DTRhoV_),ncomp,id);
397 for (auto& val : res) val = -DTRhoV_[i_it] / ( RhoV_[i_it] * RhoV_[i_it] );
398}
399
400inline void DPIRhoV_vec(SpanD T, const SpanD P, SpanD res, int ncomp, int id)
401{
402 VectorD RhoV_((int)res.size()), DPRhoV_((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
403 RhoV(T,P,SpanD(RhoV_),ncomp,id);
404 DPRhoV(T,P,SpanD(DPRhoV_),ncomp,id);
405 for (auto& val : res) val = -DPRhoV_[i_it] / ( RhoV_[i_it] * RhoV_[i_it] );
406}
407
408inline void BETAV(SpanD T, const SpanD P, SpanD res, int ncomp, int id)
409{
410 assert ((int)T.size() == ncomp * (int)res.size() && (int)T.size() == ncomp * (int)P.size());
411 VectorD IRhoV_((int)res.size()), DTIRhoV_((int)res.size()); // on peut utiliser res pour l'un de deux mais bon ...
412 IRhoV_vec(T,P,SpanD(IRhoV_),ncomp,id);
413 DTIRhoV_vec(T,P,SpanD(DTIRhoV_),ncomp,id);
414 for (auto& val : res) val = DTIRhoV_[i_it] / IRhoV_[i_it];
415}
416
417/* viscosite de la vapeur */
418inline void MuV(const SpanD T, SpanD res, int ncomp, int id)
419{
420 assert ((int)T.size() == ncomp * (int)res.size());
421 for (auto& val : res) val = 2.5e-6 + 1.5e-8 * (T[i_it * ncomp + id] + 273.15);
422}
423
424/* conductivite de la vapeur */
425inline void LambdaV(const SpanD T, SpanD res, int ncomp, int id)
426{
427 assert ((int)T.size() == ncomp * (int)res.size());
428 for (auto& val : res) val = 0.045;
429}
430
431#undef i_it
432#endif /* lois_sodium_span_included */
433
434//// TODO : FIXME : a recrire en span ... pas utilise pour le moment
435///* Lois physiques du sodium issues de fits sur DEN/DANS/DM2S/STMF/LMES/RT/12-018/A */
436//#define Tk ((T) + 273.15)
437//const double Ksil_d = 1e-9;
438///* inverses par methode de Newton */
439//inline double TLh(double h, double T0, double P)
440//{
441// double T = T0, sec;
442// int i;
443// for (i = 0; i < 100 && std::abs( sec = HL(T, P) - h ) > 1e-8; i++)
444// T -= sec / DTHL(T, P);
445// return i < 100 ? T : -1e10;
446//}
447///* energie volumique du liquide */
448//inline double EL(double T, double P) { return HL(T, P) * RhoL(T, P); }
449///* derivees */
450//inline double DTEL(double T, double P) { return DTHL(T, P) * RhoL(T, P) + HL(T, P) * DTRhoL(T, P); }
451//inline double DPEL(double T, double P) { return DPHL(T, P) * RhoL(T, P) + HL(T, P) * DPRhoL(T, P); }
452//inline double TLe(double e, double T0, double P)
453//{
454// double T = T0, sec;
455// int i;
456// for (i = 0; i < 100 && std::abs( sec = EL(T, P) - e ) > 1e-3; i++)
457// T -= sec / DTEL(T, P);
458// return i < 100 ? T : -1e10;
459//}
460///* Nusselt du liquide -> Skupinski */
461//inline double NuL(double T, double P, double w, double Dh)
462//{
463// double Pe = RhoL(T, P) * std::abs(w) * Dh * DTHL(T, P) / LambdaL(T);
464// return 4.82 + 0.0185 * pow(Pe, 0.827);
465//}
466//
467///* energie volumique de la vapeur */
468//inline double EV( double T, double P) { return RhoV(T, P) * HV(T, P); }
469///* ses derivees */
470//inline double DTEV(double T, double P) { return DTRhoV(T, P) * HV(T, P) + RhoV(T, P) * DTHV(T, P); }
471//inline double DPEV(double T, double P) { return DPRhoV(T, P) * HV(T, P) + RhoV(T, P) * DPHV(T, P); }
472//
473///* Nusselt de la vapeur -> Dittus/ Boetler */
474//inline double NuV(double T, double P, double w, double Dh)
475//{
476// double Re = RhoV(T, P) * std::abs(w) * Dh / MuV(T), Pr = MuV(T) * DTHV(T, P) / LambdaV(T);
477// return 0.023 * pow(Re, 0.8) * pow(Pr, 0.4);
478//}
479//
480//inline double Rho( double T, double x, double P ) { return 1 / ( (1-x)/RhoL(T, P) + x / RhoV(T, P) ); }