TrioCFD 1.9.8
TrioCFD documentation
Loading...
Searching...
No Matches
Fluide_R12_c1_liquide.cpp
1/****************************************************************************
2* Copyright (c) 2021, 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#include <Fluide_R12_c1_liquide.h>
17
18Implemente_instanciable(Fluide_R12_c1_liquide, "Fluide_R12_c1_liquide", Fluide_reel_base);
19
20Sortie& Fluide_R12_c1_liquide::printOn(Sortie& os) const { return os; }
21
23{
24#if HAVE_LIBC3
25 Param param(que_suis_je());
26 param.ajouter("like_eos",&like_eos_);
27 param.lire_avec_accolades_depuis(is);
28#else
29 Process::exit(que_suis_je() + " : this binary was not compiled with C3 water laws!");
30#endif
31 return is;
32}
33
34#define ind std::distance(res.begin(), &val)
35
36void Fluide_R12_c1_liquide::rho_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
37{
38 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
39#if HAVE_LIBC3
40 if (like_eos_)
41 for (auto& val : res)
42 {
43 // On fait comme dans EOS : des trucs chelous
44 // 1. Enthalpie liquide
45 int un = 1;
46 double Hl;
47 int ill, ivstat, ierrth;
48 F77NAME(FPTHLR12)(&un, &T[ind * ncomp + id], &Hl, &ill, &ivstat, &ierrth);
49
50 // 2. Saturation
51 double tsat, dtsp1 ;
52 F77NAME(FPSATR12)(&un, &P[ind], &tsat, &dtsp1, &ill, &ivstat, &ierrth);
53
54 // 3. Masse volumique a saturation
55 double DP_Tl = 0.; // Renvoye a 0. dans FCPLR12, on gange un appel
56 double rhog_sat, rhol_sat, DP_rhoL_sat, DP_rhoG_sat;
57 F77NAME(FROVLR12)(&un, &P[ind], &tsat, &tsat , &rhog_sat, &rhol_sat , &DP_Tl, &DP_rhoL_sat, &DP_rhoG_sat, &ill, &ivstat, &ierrth);
58
59 // 4. Enthalpie a saturation
60 double hvsp, hlsp, dhvsp1, dhlsp1;
61 F77NAME(FHSATR12)(&un, &tsat, &rhog_sat, &rhol_sat, &P[ind], &hvsp, &hlsp, &dtsp1, &DP_rhoG_sat, &DP_rhoL_sat, &dhvsp1, &dhlsp1);
62
63 // 5. Inversion des inconnues, sort la masse volumique...
64 double Tl, rhol;
65 F77NAME(FTLR12)(&un, &P[ind], &tsat, &rhog_sat, &hvsp , &Hl, &Tl, &rhol, &ill, &ivstat, &ierrth);
66
67 val = rhol; // NB : only a function of T
68 }
69 else
70 for (auto& val : res)
71 {
72 int un = 1;
73 int ill, ivstat, ierrth;
74
75 // 1. Saturation
76 double tsat, dtsp1 ;
77 F77NAME(FPSATR12)(&un, &P[ind], &tsat, &dtsp1, &ill, &ivstat, &ierrth);
78
79 // 2. Allez hop
80 double rhog, rhol, DP_rhoL, DP_rhoG;
81 double DP_Tl = 0.; // Renvoye a 0. dans FCPLR12, on gagne un appel
82 F77NAME(FROVLR12)(&un, &P[ind], &tsat, &T[ind * ncomp + id] , &rhog, &rhol , &DP_Tl, &DP_rhoL, &DP_rhoG, &ill, &ivstat, &ierrth);
83
84 val = rhol; // NB : only a function of T
85 }
86
87#else
88 for (auto& val : res) val = 0;
89#endif
90}
91
92void Fluide_R12_c1_liquide::dP_rho_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
93{
94 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
95#if HAVE_LIBC3
96 if (like_eos_)
97 for (auto& val : res)
98 {
99 // On fait comme dans EOS : des trucs chelous
100 // 1. Enthalpie liquide
101 int un = 1;
102 double Hl;
103 int ill, ivstat, ierrth;
104 F77NAME(FPTHLR12)(&un, &T[ind * ncomp + id], &Hl, &ill, &ivstat, &ierrth);
105
106 // 2. Saturation
107 double tsat, dtsp1 ;
108 F77NAME(FPSATR12)(&un, &P[ind], &tsat, &dtsp1, &ill, &ivstat, &ierrth);
109
110 // 3. Masse volumique a saturation
111 double DP_Tl = 0.; // Renvoye a 0. dans FCPLR12, on gange un appel
112 double rhog_sat, rhol_sat, DP_rhoL_sat, DP_rhoG_sat;
113 F77NAME(FROVLR12)(&un, &P[ind], &tsat, &tsat , &rhog_sat, &rhol_sat , &DP_Tl, &DP_rhoL_sat, &DP_rhoG_sat, &ill, &ivstat, &ierrth);
114
115 // 4. Enthalpie a saturation
116 double hvsp, hlsp, dhvsp1, dhlsp1;
117 F77NAME(FHSATR12)(&un, &tsat, &rhog_sat, &rhol_sat, &P[ind], &hvsp, &hlsp, &dtsp1, &DP_rhoG_sat, &DP_rhoL_sat, &dhvsp1, &dhlsp1);
118
119 // 5. Inversion des inconnues...
120 double Tl, rhol;
121 F77NAME(FTLR12)(&un, &P[ind], &tsat, &rhog_sat, &hvsp , &Hl, &Tl, &rhol, &ill, &ivstat, &ierrth);
122
123 // 6. Va chercher les derivees
124 double cpl, dP_Tl, dHl_Tl, dP_rhol, dHl_rhol, dP_cpl, dHl_cpl;
125 F77NAME(FCPLR12)(&un, &Tl, &cpl, &dP_Tl, &dHl_Tl, &dP_rhol, &dHl_rhol, &dP_cpl, &dHl_cpl, &ill, &ivstat, &ierrth);
126
127 val = dP_rhol - dHl_rhol* dP_Tl/dHl_Tl; // NB : this is derivative with T fixed, not with H fixed
128 }
129 else
130 for (auto& val : res)
131 {
132 int un = 1;
133 int ill, ivstat, ierrth;
134
135 // 2. Direct !!
136 double cpl, dP_Tl, dHl_Tl, dP_rhol, dHl_rhol, dP_cpl, dHl_cpl;
137 F77NAME(FCPLR12)(&un, &T[ind * ncomp + id], &cpl, &dP_Tl, &dHl_Tl, &dP_rhol, &dHl_rhol, &dP_cpl, &dHl_cpl, &ill, &ivstat, &ierrth);
138
139 }
140
141#else
142 for (auto& val : res) val = 0;
143#endif
144}
145
146void Fluide_R12_c1_liquide::dT_rho_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
147{
148 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
149#if HAVE_LIBC3
150 /* calcul a saturation */
151 if (like_eos_)
152 for (auto& val : res)
153 {
154 // On fait comme dans EOS : des trucs chelous
155 // 1. Enthalpie liquide
156 int un = 1;
157 double Hl;
158 int ill, ivstat, ierrth;
159 F77NAME(FPTHLR12)(&un, &T[ind * ncomp + id], &Hl, &ill, &ivstat, &ierrth);
160
161 // 2. Saturation
162 double tsat, dtsp1 ;
163 F77NAME(FPSATR12)(&un, &P[ind], &tsat, &dtsp1, &ill, &ivstat, &ierrth);
164
165 // 3. Masse volumique a saturation
166 double DP_Tl = 0.; // Renvoye a 0. dans FCPLR12, on gange un appel
167 double rhog_sat, rhol_sat, DP_rhoL_sat, DP_rhoG_sat;
168 F77NAME(FROVLR12)(&un, &P[ind], &tsat, &tsat , &rhog_sat, &rhol_sat , &DP_Tl, &DP_rhoL_sat, &DP_rhoG_sat, &ill, &ivstat, &ierrth);
169
170 // 4. Enthalpie a saturation
171 double hvsp, hlsp, dhvsp1, dhlsp1;
172 F77NAME(FHSATR12)(&un, &tsat, &rhog_sat, &rhol_sat, &P[ind], &hvsp, &hlsp, &dtsp1, &DP_rhoG_sat, &DP_rhoL_sat, &dhvsp1, &dhlsp1);
173
174 // 5. Inversion des inconnues...
175 double Tl, rhol;
176 F77NAME(FTLR12)(&un, &P[ind], &tsat, &rhog_sat, &hvsp , &Hl, &Tl, &rhol, &ill, &ivstat, &ierrth);
177
178 // 6. Va chercher les derivees
179 double cpl, dP_Tl, dHl_Tl, dP_rhol, dHl_rhol, dP_cpl, dHl_cpl;
180 F77NAME(FCPLR12)(&un, &Tl, &cpl, &dP_Tl, &dHl_Tl, &dP_rhol, &dHl_rhol, &dP_cpl, &dHl_cpl, &ill, &ivstat, &ierrth);
181
182 val = dHl_rhol/dHl_Tl;
183 }
184 else for (auto& val : res)
185 {
186 // 1. Si facile
187 int un = 1;
188 int ill, ivstat, ierrth;
189 double cpl, dP_Tl, dHl_Tl, dP_rhol, dHl_rhol, dP_cpl, dHl_cpl;
190 F77NAME(FCPLR12)(&un, &T[ind * ncomp + id], &cpl, &dP_Tl, &dHl_Tl, &dP_rhol, &dHl_rhol, &dP_cpl, &dHl_cpl, &ill, &ivstat, &ierrth);
191
192 val = dHl_rhol/dHl_Tl;
193 }
194#else
195 for (auto& val : res) val = 0;
196#endif
197}
198
199void Fluide_R12_c1_liquide::h_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
200{
201 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
202#if HAVE_LIBC3
203 /* calcul a saturation */
204 for (auto& val : res)
205 {
206 int un = 1;
207 double hl; //sortie
208 int ill, ivstat, ierrth;
209 F77NAME(FPTHLR12)(&un, &T[ind * ncomp + id], &hl, &ill, &ivstat, &ierrth);
210 val = hl;
211 }
212#else
213 for (auto& val : res) val = 0;
214#endif
215}
216
217void Fluide_R12_c1_liquide::dP_h_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
218{
219 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
220#if HAVE_LIBC3
221 if (like_eos_)
222 for (auto& val : res)
223 {
224 // On fait comme dans EOS : des trucs chelous
225 // 1. Enthalpie liquide
226 int un = 1;
227 double Hl;
228 int ill, ivstat, ierrth;
229 F77NAME(FPTHLR12)(&un, &T[ind * ncomp + id], &Hl, &ill, &ivstat, &ierrth);
230
231 // 2. Saturation
232 double tsat, dtsp1 ;
233 F77NAME(FPSATR12)(&un, &P[ind], &tsat, &dtsp1, &ill, &ivstat, &ierrth);
234
235 // 3. Masse volumique a saturation
236 double DP_Tl = 0.; // Renvoye a 0. dans FCPLR12, on gange un appel
237 double rhog_sat, rhol_sat, DP_rhoL_sat, DP_rhoG_sat;
238 F77NAME(FROVLR12)(&un, &P[ind], &tsat, &tsat , &rhog_sat, &rhol_sat , &DP_Tl, &DP_rhoL_sat, &DP_rhoG_sat, &ill, &ivstat, &ierrth);
239
240 // 4. Enthalpie a saturation
241 double hvsp, hlsp, dhvsp1, dhlsp1;
242 F77NAME(FHSATR12)(&un, &tsat, &rhog_sat, &rhol_sat, &P[ind], &hvsp, &hlsp, &dtsp1, &DP_rhoG_sat, &DP_rhoL_sat, &dhvsp1, &dhlsp1);
243
244 // 5. Inversion des inconnues...
245 double Tl, rhol;
246 F77NAME(FTLR12)(&un, &P[ind], &tsat, &rhog_sat, &hvsp , &Hl, &Tl, &rhol, &ill, &ivstat, &ierrth);
247
248 // 6. Va chercher les derivees
249 double cpl, dP_Tl, dHl_Tl, dP_rhol, dHl_rhol, dP_cpl, dHl_cpl;
250 F77NAME(FCPLR12)(&un, &Tl, &cpl, &dP_Tl, &dHl_Tl, &dP_rhol, &dHl_rhol, &dP_cpl, &dHl_cpl, &ill, &ivstat, &ierrth);
251
252 val = -dP_Tl/dHl_Tl;
253 }
254 else for (auto& val : res)
255 {
256 int un = 1;
257 int ill, ivstat, ierrth;
258 // 1. Si beau
259 double cpl, dP_Tl, dHl_Tl, dP_rhol, dHl_rhol, dP_cpl, dHl_cpl;
260 F77NAME(FCPLR12)(&un, &T[ind * ncomp + id], &cpl, &dP_Tl, &dHl_Tl, &dP_rhol, &dHl_rhol, &dP_cpl, &dHl_cpl, &ill, &ivstat, &ierrth);
261
262 val = -dP_Tl/dHl_Tl;
263 }
264#else
265 for (auto& val : res) val = 0;
266#endif
267}
268
269void Fluide_R12_c1_liquide::dT_h_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
270{
271 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
272#if HAVE_LIBC3
273 /* calcul a saturation */
274 if (like_eos_)
275 for (auto& val : res)
276 {
277 // On fait comme dans EOS : des trucs chelous
278 // 1. Enthalpie liquide
279 int un = 1;
280 double Hl;
281 int ill, ivstat, ierrth;
282 F77NAME(FPTHLR12)(&un, &T[ind * ncomp + id], &Hl, &ill, &ivstat, &ierrth);
283
284 // 2. Saturation
285 double tsat, dtsp1 ;
286 F77NAME(FPSATR12)(&un, &P[ind], &tsat, &dtsp1, &ill, &ivstat, &ierrth);
287
288 // 3. Masse volumique a saturation
289 double DP_Tl = 0.; // Renvoye a 0. dans FCPLR12, on gange un appel
290 double rhog_sat, rhol_sat, DP_rhoL_sat, DP_rhoG_sat;
291 F77NAME(FROVLR12)(&un, &P[ind], &tsat, &tsat , &rhog_sat, &rhol_sat , &DP_Tl, &DP_rhoL_sat, &DP_rhoG_sat, &ill, &ivstat, &ierrth);
292
293 // 4. Enthalpie a saturation
294 double hvsp, hlsp, dhvsp1, dhlsp1;
295 F77NAME(FHSATR12)(&un, &tsat, &rhog_sat, &rhol_sat, &P[ind], &hvsp, &hlsp, &dtsp1, &DP_rhoG_sat, &DP_rhoL_sat, &dhvsp1, &dhlsp1);
296
297 // 5. Inversion des inconnues...
298 double Tl, rhol;
299 F77NAME(FTLR12)(&un, &P[ind], &tsat, &rhog_sat, &hvsp , &Hl, &Tl, &rhol, &ill, &ivstat, &ierrth);
300
301 // 6. Va chercher les derivees
302 double cpl, dP_Tl, dHl_Tl, dP_rhol, dHl_rhol, dP_cpl, dHl_cpl;
303 F77NAME(FCPLR12)(&un, &Tl, &cpl, &dP_Tl, &dHl_Tl, &dP_rhol, &dHl_rhol, &dP_cpl, &dHl_cpl, &ill, &ivstat, &ierrth);
304
305 val = 1./dHl_Tl; // =cpl normalement... mais ici non !
306
307
308 // 7. On teste la version plus rapide
309 double cpl2, dP_Tl2, dHl_Tl2, dP_rhol2, dHl_rhol2, dP_cpl2, dHl_cpl2;
310 F77NAME(FCPLR12)(&un, &T[ind * ncomp + id], &cpl2, &dP_Tl2, &dHl_Tl2, &dP_rhol2, &dHl_rhol2, &dP_cpl2, &dHl_cpl2, &ill, &ivstat, &ierrth);
311 if ((abs(1./dHl_Tl)>1.e-10) && (abs((1./dHl_Tl-1./dHl_Tl2)/(1./dHl_Tl)) > 1.e-4 )) Cout << "dT_h2 : " << 1./dHl_Tl2 << " dT_h : " << 1./dHl_Tl << finl ;
312 }
313 else for (auto& val : res)
314 {
315 int un = 1;
316 int ill, ivstat, ierrth;
317 double cpl, dP_Tl, dHl_Tl, dP_rhol, dHl_rhol, dP_cpl, dHl_cpl;
318 F77NAME(FCPLR12)(&un, &T[ind * ncomp + id], &cpl, &dP_Tl, &dHl_Tl, &dP_rhol, &dHl_rhol, &dP_cpl, &dHl_cpl, &ill, &ivstat, &ierrth);
319 val = 1./dHl_Tl; // =cpl normalement... mais ici non !
320 }
321#else
322 for (auto& val : res) val = 0;
323#endif
324}
325
326void Fluide_R12_c1_liquide::cp_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
327{
328 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
329#if HAVE_LIBC3
330 /* calcul a saturation */
331 for (auto& val : res)
332 {
333 int un = 1;
334 double cpl, dP_Tl, DHl_Tl, DP_Rhol, DHl_Rhol, DP_cpl, DHl_cpl ; //sorties
335 int ill, ivstat, ierrth;
336 F77NAME(FCPLR12)(&un, &T[ind * ncomp + id], &cpl, &dP_Tl, &DHl_Tl, &DP_Rhol, &DHl_Rhol, &DP_cpl, &DHl_cpl, &ill, &ivstat, &ierrth);
337 val = cpl;
338 }
339#else
340 for (auto& val : res) val = 0;
341#endif
342}
343
344void Fluide_R12_c1_liquide::beta_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
345{
346 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
347 VectorD dT_rho___((int )res.size()), rho___((int )res.size());
348 dT_rho_(T,P,SpanD(dT_rho___),ncomp,id);
349 rho_(T,P,SpanD(rho___),ncomp,id);
350 for (auto& val : res) val = dT_rho___[ind] / rho___[ind];
351}
352
353void Fluide_R12_c1_liquide::mu_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
354{
355 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
356#if HAVE_LIBC3
357 /* calcul a saturation */
358 for (auto& val : res)
359 {
360 // Fit we did manually using NIST data
361 double a1 = 4.51272790e-09;
362 double a2 = -2.14503254e-06;
363 double a3 = 2.44480131e-04 ;
364
365 val = a3 + a2*T[ind * ncomp + id] + a1*T[ind * ncomp + id]*T[ind * ncomp + id];
366 }
367#else
368 for (auto& val : res) val = 0;
369#endif
370}
371
372void Fluide_R12_c1_liquide::lambda_(const SpanD T, const SpanD P, SpanD res, int ncomp, int id) const
373{
374 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )res.size());
375#if HAVE_LIBC3
376 /* calcul a saturation */
377 for (auto& val : res)
378 {
379 int un = 1;
380 double dtl1 = 0., dtl2 = 0.; // We don't need the derivative of lambda
381 double lambdal, dlambl1, dlambl2;
382 F77NAME(FCONLR12)(&un, &T[ind * ncomp + id], &dtl1, &dtl2, &lambdal, &dlambl1, &dlambl2);
383 val = lambdal;
384 }
385#else
386 for (auto& val : res) val = 0;
387#endif
388}
389
390#undef ind
Class defining operators and methods for all reading operation in an input flow (file,...
Definition Entree.h:42
void h_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void dT_rho_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void rho_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void mu_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void dP_rho_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void dT_h_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void dP_h_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void beta_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void lambda_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
void cp_(const SpanD T, const SpanD P, SpanD res, int ncomp=1, int id=0) const override
Classe Fluide_reel_base Cette classe represente un fluide reel ainsi que.
const Nom & que_suis_je() const
renvoie la chaine identifiant la classe.
Definition Objet_U.cpp:104
virtual Entree & readOn(Entree &)
Lecture d'un Objet_U sur un flot d'entree Methode a surcharger.
Definition Objet_U.cpp:293
virtual Sortie & printOn(Sortie &) const
Ecriture de l'objet sur un flot de sortie Methode a surcharger.
Definition Objet_U.cpp:282
static void exit(int exit_code=-1)
Routine de sortie de TRUST dans une region Kokkos.
Definition Process.cpp:455
Classe de base des flux de sortie.
Definition Sortie.h:52