TrioCFD 1.9.8
TrioCFD documentation
Loading...
Searching...
No Matches
CoolProp_to_TRUST_generique.cpp
1/****************************************************************************
2* Copyright (c) 2024, 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 <CoolProp_to_TRUST_generique.h>
17
18using namespace CoolProp;
19
20// iterator index !
21#define i_it std::distance(TT.begin(), &val)
22#define bi_it std::distance(bTT.begin(), &bval)
23
24void CoolProp_to_TRUST_generique::set_fluide_generique(const char *const model_name, const char *const fluid_name)
25{
26#ifdef HAS_COOLPROP
27 fluide = std::shared_ptr<AbstractState>(AbstractState::factory(std::string(model_name), std::string(fluid_name)));
28#else
29 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
30 throw;
31#endif
32}
33
34int CoolProp_to_TRUST_generique::tppi_get_single_property_T_(Loi_en_T enum_prop, const SpanD P, const SpanD T, SpanD R, int ncomp, int ind) const
35{
36#ifdef HAS_COOLPROP
37 assert((int )T.size() == ncomp * (int )P.size() && (int )T.size() == ncomp * (int )R.size());
38 if (ncomp == 1) return tppi_get_single_property_T__(enum_prop, P, T, R);
39 else /* attention stride */
40 {
41 VectorD temp_((int)P.size());
42 SpanD TT(temp_);
43 for (auto& val : TT) val = T[i_it * ncomp + ind];
44 return tppi_get_single_property_T__(enum_prop, P, TT, R);
45 }
46#else
47 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
48 throw;
49#endif
50}
51
52int CoolProp_to_TRUST_generique::tppi_get_single_property_T__(Loi_en_T enum_prop, const SpanD P, const SpanD T, SpanD R) const
53{
54#ifdef HAS_COOLPROP
55 const int sz = (int )R.size();
56 // derivees qui manquent ...
57 if (enum_prop == Loi_en_T::MU_DP || enum_prop == Loi_en_T::LAMBDA_DP || enum_prop == Loi_en_T::SIGMA_DP)
58 return FD_derivative_pT__(enum_prop,P, T, R);
59
60 if (enum_prop == Loi_en_T::MU_DT || enum_prop == Loi_en_T::LAMBDA_DT || enum_prop == Loi_en_T::SIGMA_DT)
61 return FD_derivative_pT__(enum_prop,P, T, R, false /* wrt T */);
62
63 for (int i = 0; i < sz; i++)
64 {
65 fluide->update(CoolProp::PT_INPUTS, P[i], T[i]); // SI units
66
67 if (enum_prop == Loi_en_T::RHO) R[i] = fluide->rhomass();
68 if (enum_prop == Loi_en_T::RHO_DP) R[i] = fluide->first_partial_deriv(CoolProp::iDmass, CoolProp::iP, CoolProp::iT);
69 if (enum_prop == Loi_en_T::RHO_DT) R[i] = fluide->first_partial_deriv(CoolProp::iDmass, CoolProp::iT, CoolProp::iP);
70
71 if (enum_prop == Loi_en_T::H) R[i] = fluide->hmass();
72 if (enum_prop == Loi_en_T::H_DP) R[i] = fluide->first_partial_deriv(CoolProp::iHmass, CoolProp::iP, CoolProp::iT);
73 if (enum_prop == Loi_en_T::H_DT) R[i] = fluide->first_partial_deriv(CoolProp::iHmass, CoolProp::iT, CoolProp::iP);
74
75 if (enum_prop == Loi_en_T::CP) R[i] = fluide->cpmass();
76 if (enum_prop == Loi_en_T::CP_DP) R[i] = fluide->first_partial_deriv(CoolProp::iCpmass, CoolProp::iP, CoolProp::iT);
77 if (enum_prop == Loi_en_T::CP_DT) R[i] = fluide->first_partial_deriv(CoolProp::iCpmass, CoolProp::iT, CoolProp::iP);
78
79 if (enum_prop == Loi_en_T::MU) R[i] = fluide->viscosity();
80 if (enum_prop == Loi_en_T::LAMBDA) R[i] = fluide->conductivity();
81 if (enum_prop == Loi_en_T::BETA) R[i] = fluide->isobaric_expansion_coefficient();
82 if (enum_prop == Loi_en_T::SIGMA) R[i] = fluide->surface_tension();
83 }
84 return 0; // FIXME : on suppose que tout OK
85#else
86 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
87 throw;
88#endif
89}
90
91int CoolProp_to_TRUST_generique::FD_derivative_pT__(Loi_en_T enum_prop, const SpanD P, const SpanD T, SpanD R, bool wrt_p) const
92{
93#ifdef HAS_COOLPROP
94 const int sz = (int )R.size();
95 for (int i = 0; i < sz; i++)
96 {
97 double plus_ = EPS, minus_ = EPS;
98
99 wrt_p ? fluide->update(CoolProp::PT_INPUTS, P[i] * (1. + EPS), T[i]) : fluide->update(CoolProp::PT_INPUTS, P[i], T[i] * (1. + EPS));
100
101 if (enum_prop == Loi_en_T::MU_DP || enum_prop == Loi_en_T::MU_DT) plus_ = fluide->viscosity();
102 if (enum_prop == Loi_en_T::LAMBDA_DP || enum_prop == Loi_en_T::LAMBDA_DT) plus_ = fluide->conductivity();
103 if (enum_prop == Loi_en_T::SIGMA_DP || enum_prop == Loi_en_T::SIGMA_DT) plus_ = fluide->surface_tension();
104
105 wrt_p ? fluide->update(CoolProp::PT_INPUTS, P[i] * (1. - EPS), T[i]) : fluide->update(CoolProp::PT_INPUTS, P[i], T[i] * (1. - EPS));
106
107 if (enum_prop == Loi_en_T::MU_DP || enum_prop == Loi_en_T::MU_DT) minus_ = fluide->viscosity();
108 if (enum_prop == Loi_en_T::LAMBDA_DP || enum_prop == Loi_en_T::LAMBDA_DT) minus_ = fluide->conductivity();
109 if (enum_prop == Loi_en_T::SIGMA_DP || enum_prop == Loi_en_T::SIGMA_DT) minus_ = fluide->surface_tension();
110
111 R[i] = wrt_p ? ((plus_ - minus_) / ( 2 * EPS * P[i])) : ((plus_ - minus_) / ( 2 * EPS * T[i]));
112 }
113 return 0; // FIXME : on suppose que tout OK
114#else
115 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
116 throw;
117#endif
118}
119
120int CoolProp_to_TRUST_generique::tppi_get_single_property_h_(Loi_en_h enum_prop, const SpanD P, const SpanD H, SpanD R, int ncomp, int ind) const
121{
122#ifdef HAS_COOLPROP
123 assert((int )H.size() == ncomp * (int )P.size() && (int )H.size() == ncomp * (int )R.size());
124 if (ncomp == 1) return tppi_get_single_property_h__(enum_prop, P, H, R);
125 else /* attention stride */
126 {
127 VectorD temp_((int)P.size());
128 SpanD TT(temp_);
129 for (auto& val : TT) val = H[i_it * ncomp + ind];
130 return tppi_get_single_property_h__(enum_prop, P, TT, R);
131 }
132 return 0; // FIXME : on suppose que tout OK
133#else
134 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
135 throw;
136#endif
137}
138
139int CoolProp_to_TRUST_generique::tppi_get_single_property_h__(Loi_en_h enum_prop, const SpanD P, const SpanD H, SpanD R) const
140{
141#ifdef HAS_COOLPROP
142 const int sz = (int )R.size();
143 // derivees qui manquent ...
144 if (enum_prop == Loi_en_h::MU_DP || enum_prop == Loi_en_h::LAMBDA_DP || enum_prop == Loi_en_h::SIGMA_DP)
145 return FD_derivative_ph__(enum_prop,P, H, R);
146
147 if (enum_prop == Loi_en_h::MU_DH || enum_prop == Loi_en_h::LAMBDA_DH || enum_prop == Loi_en_h::SIGMA_DH)
148 return FD_derivative_ph__(enum_prop,P, H, R, false /* wrt H */);
149
150 for (int i = 0; i < sz; i++)
151 {
152 fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i]); // SI units
153
154 if (enum_prop == Loi_en_h::RHO) R[i] = fluide->rhomass();
155 if (enum_prop == Loi_en_h::RHO_DP) R[i] = fluide->first_partial_deriv(CoolProp::iDmass, CoolProp::iP, CoolProp::iHmass);
156 if (enum_prop == Loi_en_h::RHO_DH) R[i] = fluide->first_partial_deriv(CoolProp::iDmass, CoolProp::iHmass, CoolProp::iP);
157
158 if (enum_prop == Loi_en_h::T) R[i] = fluide->T();
159 if (enum_prop == Loi_en_h::T_DP) R[i] = fluide->first_partial_deriv(CoolProp::iT, CoolProp::iP, CoolProp::iHmass);
160 if (enum_prop == Loi_en_h::T_DH) R[i] = fluide->first_partial_deriv(CoolProp::iT, CoolProp::iHmass, CoolProp::iP);
161
162 if (enum_prop == Loi_en_h::CP) R[i] = fluide->cpmass();
163 if (enum_prop == Loi_en_h::CP_DP) R[i] = fluide->first_partial_deriv(CoolProp::iCpmass, CoolProp::iP, CoolProp::iHmass);
164 if (enum_prop == Loi_en_h::CP_DH) R[i] = fluide->first_partial_deriv(CoolProp::iCpmass, CoolProp::iHmass, CoolProp::iP);
165
166 if (enum_prop == Loi_en_h::MU) R[i] = fluide->viscosity();
167 if (enum_prop == Loi_en_h::LAMBDA) R[i] = fluide->conductivity();
168 if (enum_prop == Loi_en_h::BETA) R[i] = fluide->isobaric_expansion_coefficient();
169 if (enum_prop == Loi_en_h::SIGMA) R[i] = fluide->surface_tension();
170 }
171 return 0; // FIXME : on suppose que tout OK
172#else
173 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
174 throw;
175#endif
176}
177
178int CoolProp_to_TRUST_generique::FD_derivative_ph__(Loi_en_h enum_prop, const SpanD P, const SpanD H, SpanD R, bool wrt_p) const
179{
180#ifdef HAS_COOLPROP
181 const int sz = (int )R.size();
182 for (int i = 0; i < sz; i++)
183 {
184 double plus_ = EPS, minus_ = EPS;
185
186 wrt_p ? fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i] * (1. + EPS)) : fluide->update(CoolProp::HmassP_INPUTS, H[i] * (1. + EPS), P[i]);
187
188 if (enum_prop == Loi_en_h::MU_DP || enum_prop == Loi_en_h::MU_DH) plus_ = fluide->viscosity();
189 if (enum_prop == Loi_en_h::LAMBDA_DP || enum_prop == Loi_en_h::LAMBDA_DH) plus_ = fluide->conductivity();
190 if (enum_prop == Loi_en_h::SIGMA_DP || enum_prop == Loi_en_h::SIGMA_DH) plus_ = fluide->surface_tension();
191
192 wrt_p ? fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i] * (1. - EPS)) : fluide->update(CoolProp::HmassP_INPUTS, H[i] * (1. - EPS), P[i]);
193
194 if (enum_prop == Loi_en_h::MU_DP || enum_prop == Loi_en_h::MU_DH) minus_ = fluide->viscosity();
195 if (enum_prop == Loi_en_h::LAMBDA_DP || enum_prop == Loi_en_h::LAMBDA_DH) minus_ = fluide->conductivity();
196 if (enum_prop == Loi_en_h::SIGMA_DP || enum_prop == Loi_en_h::SIGMA_DH) minus_ = fluide->surface_tension();
197
198 R[i] = wrt_p ? ((plus_ - minus_) / ( 2 * EPS * P[i])) : ((plus_ - minus_) / ( 2 * EPS * H[i]));
199 }
200 return 0; // FIXME : on suppose que tout OK
201#else
202 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
203 throw;
204#endif
205}
206
207int CoolProp_to_TRUST_generique::tppi_get_all_properties_T__(const MSpanD input, MLoiSpanD prop) const
208{
209#ifdef HAS_COOLPROP
210 const SpanD T = input.at("temperature"), P = input.at("pressure");
211 const int sz = (int ) P.size();
212 for (int i = 0; i < sz; i++)
213 {
214 fluide->update(CoolProp::PT_INPUTS, P[i], T[i]); // SI units
215 for (auto &itr : prop)
216 {
217 Loi_en_T prop_ = itr.first;
218 SpanD span_ = itr.second;
219
220 if (prop_ == Loi_en_T::RHO) span_[i] = fluide->rhomass();
221 else if (prop_ == Loi_en_T::RHO_DP) span_[i] = fluide->first_partial_deriv(CoolProp::iDmass, CoolProp::iP, CoolProp::iT);
222 else if (prop_ == Loi_en_T::RHO_DT) span_[i] = fluide->first_partial_deriv(CoolProp::iDmass, CoolProp::iT, CoolProp::iP);
223 else if (prop_ == Loi_en_T::H) span_[i] = fluide->hmass();
224 else if (prop_ == Loi_en_T::H_DP) span_[i] = fluide->first_partial_deriv(CoolProp::iHmass, CoolProp::iP, CoolProp::iT);
225 else if (prop_ == Loi_en_T::H_DT) span_[i] = fluide->first_partial_deriv(CoolProp::iHmass, CoolProp::iT, CoolProp::iP);
226 else if (prop_ == Loi_en_T::CP) span_[i] = fluide->cpmass();
227 else if (prop_ == Loi_en_T::CP_DP) span_[i] = fluide->first_partial_deriv(CoolProp::iCpmass, CoolProp::iP, CoolProp::iT);
228 else if (prop_ == Loi_en_T::CP_DT) span_[i] = fluide->first_partial_deriv(CoolProp::iCpmass, CoolProp::iT, CoolProp::iP);
229 else if (prop_ == Loi_en_T::MU) span_[i] = fluide->viscosity();
230 else if (prop_ == Loi_en_T::LAMBDA) span_[i] = fluide->conductivity();
231 else if (prop_ == Loi_en_T::BETA) span_[i] = fluide->isobaric_expansion_coefficient();
232 else if (prop_ == Loi_en_T::SIGMA) span_[i] = fluide->surface_tension();
233 else Process::exit("Derivative requested is not coded !");
234 }
235 }
236 return 0; // FIXME : on suppose que tout OK
237#else
238 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
239 throw;
240#endif
241}
242
243int CoolProp_to_TRUST_generique::tppi_get_all_properties_T_IF97__(const MSpanD input, MLoiSpanD prop) const
244{
245#ifdef HAS_COOLPROP
246 const SpanD T = input.at("temperature"), P = input.at("pressure");
247 const int sz = (int ) P.size();
248
249 bool has_prop = false;
250 for (auto &itr : prop)
251 if (itr.first == Loi_en_T::RHO || itr.first == Loi_en_T::H || itr.first == Loi_en_T::CP||
252 itr.first == Loi_en_T::MU || itr.first == Loi_en_T::LAMBDA || itr.first == Loi_en_T::BETA || itr.first == Loi_en_T::SIGMA)
253 {
254 has_prop = true;
255 break;
256 }
257
258 bool has_beta = false;
259 for (auto &itr : prop)
260 if (itr.first == Loi_en_T::BETA)
261 {
262 has_beta = true;
263 break;
264 }
265
266 if (has_prop)
267 for (int i = 0; i < sz; i++)
268 {
269 fluide->update(CoolProp::PT_INPUTS, P[i], T[i]); // SI units
270 for (auto &itr : prop)
271 {
272 Loi_en_T prop_ = itr.first;
273 SpanD span_ = itr.second;
274
275 if (prop_ == Loi_en_T::RHO) span_[i] = fluide->rhomass();
276 if (prop_ == Loi_en_T::H) span_[i] = fluide->hmass();
277 if (prop_ == Loi_en_T::CP) span_[i] = fluide->cpmass();
278 if (prop_ == Loi_en_T::MU) span_[i] = fluide->viscosity();
279 if (prop_ == Loi_en_T::LAMBDA) span_[i] = fluide->conductivity();
280 if (prop_ == Loi_en_T::SIGMA) span_[i] = fluide->surface_tension();
281
282 // XXX : Elie Saikali : Attention, pas possible d'appeler beta directement car implementee en derivee dans coolprop ... on fait a la main !
283// if (prop_ == Loi_en_T::BETA) span_[i] = fluide->isobaric_expansion_coefficient();
284 if (prop_ == Loi_en_T::BETA && has_beta)
285 {
286 const double rho_ = fluide->rhomass();
287 fluide->update(CoolProp::PT_INPUTS, P[i], T[i] * (1. + EPS));
288 const double plus_rho_ = fluide->rhomass();
289 fluide->update(CoolProp::PT_INPUTS, P[i], T[i] * (1. - EPS));
290 const double minus_rho_ = fluide->rhomass();
291 const double drho_dt_ = (plus_rho_ - minus_rho_) / ( 2 * EPS * T[i]);
292 span_[i] = drho_dt_ / rho_;
293 }
294 }
295 }
296
297 // derivees qui manquent ...
298 bool has_DP_DH = false;
299 for (auto &itr : prop)
300 if (itr.first == Loi_en_T::RHO_DP || itr.first == Loi_en_T::H_DP || itr.first == Loi_en_T::CP_DP ||
301 itr.first == Loi_en_T::RHO_DT || itr.first == Loi_en_T::H_DT || itr.first == Loi_en_T::CP_DT)
302 {
303 has_DP_DH = true;
304 break;
305 }
306
307 if (has_DP_DH)
308 for (int i = 0; i < sz; i++)
309 {
310 // en P
311 double plus_rho_ = EPS, plus_H_ = EPS, plus_cp_ = EPS;
312 fluide->update(CoolProp::PT_INPUTS, P[i] * (1. + EPS), T[i]);
313 for (auto &itr : prop)
314 {
315 Loi_en_T prop_ = itr.first;
316 if (prop_ == Loi_en_T::RHO_DP) plus_rho_ = fluide->rhomass();
317 if (prop_ == Loi_en_T::H_DP) plus_H_ = fluide->hmass();
318 if (prop_ == Loi_en_T::CP_DP) plus_cp_ = fluide->cpmass();
319 }
320
321 double minus_rho_ = EPS, minus_H_ = EPS, minus_cp_ = EPS;
322 fluide->update(CoolProp::PT_INPUTS, P[i] * (1. - EPS), T[i]);
323 for (auto &itr : prop)
324 {
325 Loi_en_T prop_ = itr.first;
326 if (prop_ == Loi_en_T::RHO_DP) minus_rho_ = fluide->rhomass();
327 if (prop_ == Loi_en_T::H_DP) minus_H_ = fluide->hmass();
328 if (prop_ == Loi_en_T::CP_DP) minus_cp_ = fluide->cpmass();
329 }
330
331 for (auto &itr : prop)
332 {
333 Loi_en_T prop_ = itr.first;
334 SpanD span_ = itr.second;
335
336 if (prop_ == Loi_en_T::RHO_DP) span_[i] = (plus_rho_ - minus_rho_) / ( 2 * EPS * P[i]);
337 if (prop_ == Loi_en_T::H_DP) span_[i] = (plus_H_ - minus_H_) / ( 2 * EPS * P[i]);
338 if (prop_ == Loi_en_T::CP_DP) span_[i] = (plus_cp_ - minus_cp_) / ( 2 * EPS * P[i]);
339 }
340
341 // en T
342 plus_rho_ = EPS, plus_H_ = EPS, plus_cp_ = EPS;
343 fluide->update(CoolProp::PT_INPUTS, P[i], T[i] * (1. + EPS));
344 for (auto &itr : prop)
345 {
346 Loi_en_T prop_ = itr.first;
347 if (prop_ == Loi_en_T::RHO_DT) plus_rho_ = fluide->rhomass();
348 if (prop_ == Loi_en_T::H_DT) plus_H_ = fluide->hmass();
349 if (prop_ == Loi_en_T::CP_DT) plus_cp_ = fluide->cpmass();
350 }
351
352 minus_rho_ = EPS, minus_H_ = EPS, minus_cp_ = EPS;
353 fluide->update(CoolProp::PT_INPUTS, P[i], T[i] * (1. - EPS));
354 for (auto &itr : prop)
355 {
356 Loi_en_T prop_ = itr.first;
357 if (prop_ == Loi_en_T::RHO_DT) minus_rho_ = fluide->rhomass();
358 if (prop_ == Loi_en_T::H_DT) minus_H_ = fluide->hmass();
359 if (prop_ == Loi_en_T::CP_DT) minus_cp_ = fluide->cpmass();
360 }
361
362 for (auto &itr : prop)
363 {
364 Loi_en_T prop_ = itr.first;
365 SpanD span_ = itr.second;
366
367 if (prop_ == Loi_en_T::RHO_DT) span_[i] = (plus_rho_ - minus_rho_) / ( 2 * EPS * T[i]);
368 if (prop_ == Loi_en_T::H_DT) span_[i] = (plus_H_ - minus_H_) / ( 2 * EPS * T[i]);
369 if (prop_ == Loi_en_T::CP_DT) span_[i] = (plus_cp_ - minus_cp_) / ( 2 * EPS * T[i]);
370 }
371 }
372
373 return 0; // FIXME : on suppose que tout OK
374#else
375 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
376 throw;
377#endif
378}
379
380// methodes particulieres par application pour gagner en performance : utilisees dans Pb_Multiphase (pour le moment !)
381int CoolProp_to_TRUST_generique::tppi_get_CPMLB_pb_multiphase_pT(const MSpanD input, MLoiSpanD prop, int ncomp, int ind) const
382{
383#ifdef HAS_COOLPROP
384 const SpanD T = input.at("temperature"), P = input.at("pressure");
385
386#ifndef NDEBUG
387 assert((int )prop.size() == 4 && (int )input.size() == 2);
388 assert((int )T.size() == ncomp * (int )P.size());
389 for (auto& itr : prop) assert((int )T.size() == ncomp * (int )itr.second.size());
390#endif
391
392 Tk_(T); // XXX : ATTENTION : need Kelvin
393
394 if (ncomp == 1)
395 (fluide->backend_name() == "IF97Backend") ? tppi_get_all_properties_T_IF97__(input,prop) : tppi_get_all_properties_T__(input,prop);
396 else /* attention stride */
397 {
398 VectorD temp_((int)P.size());
399 SpanD TT(temp_);
400 for (auto& val : TT) val = T[i_it * ncomp + ind];
401 MSpanD input_ = { { "temperature", TT }, { "pressure", P }};
402
403 (fluide->backend_name() == "IF97Backend") ? tppi_get_all_properties_T_IF97__(input_,prop) : tppi_get_all_properties_T__(input_,prop);
404 }
405
406 Tc_(T); // XXX : ATTENTION : need to put back T in C
407
408 return 0; // FIXME : on suppose que tout OK
409#else
410 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
411 throw;
412#endif
413}
414
415int CoolProp_to_TRUST_generique::tppi_get_all_pb_multiphase_pT(const MSpanD input, MLoiSpanD inter, MLoiSpanD bord, int ncomp, int ind) const
416{
417#ifdef HAS_COOLPROP
418 const SpanD T = input.at("temperature"), P = input.at("pressure"), bT = input.at("bord_temperature"), bP = input.at("bord_pressure");
419
420#ifndef NDEBUG
421 assert( (int )input.size() == 4 && (int )inter.size() == 6 && (int )bord.size() == 2);
422 assert ((int )bT.size() == ncomp * (int )bP.size() && (int )T.size() == ncomp * (int )P.size());
423 for (auto& itr : inter) assert((int )T.size() == ncomp * (int )itr.second.size());
424 for (auto& itr : bord) assert((int )bT.size() == ncomp * (int )itr.second.size());
425#endif
426
427 Tk_(T), Tk_(bT); // XXX : ATTENTION : need Kelvin
428
429 if (ncomp == 1)
430 {
431 MSpanD input_int_ = { { "temperature", T }, { "pressure", P }};
432 MSpanD input_bord_ = { { "temperature", bT }, { "pressure", bP }};
433 if ((fluide->backend_name() == "IF97Backend"))
434 {
435 tppi_get_all_properties_T_IF97__(input_int_,inter);
436 tppi_get_all_properties_T_IF97__(input_bord_,bord);
437 }
438 else
439 {
440 tppi_get_all_properties_T__(input_int_,inter);
441 tppi_get_all_properties_T__(input_bord_,bord);
442 }
443 }
444 else /* attention stride */
445 {
446 VectorD temp_((int)P.size()), btemp_((int)bP.size());
447 SpanD TT(temp_), bTT(btemp_);
448 for (auto& val : TT) val = T[i_it * ncomp + ind];
449 for (auto& bval : bTT) bval = bT[bi_it * ncomp + ind];
450
451 MSpanD input_int_ = { { "temperature", TT }, { "pressure", P }};
452 MSpanD input_bord_ = { { "temperature", bTT }, { "pressure", bP }};
453 if ((fluide->backend_name() == "IF97Backend"))
454 {
455 tppi_get_all_properties_T_IF97__(input_int_,inter);
456 tppi_get_all_properties_T_IF97__(input_bord_,bord);
457 }
458 else
459 {
460 tppi_get_all_properties_T__(input_int_,inter);
461 tppi_get_all_properties_T__(input_bord_,bord);
462 }
463 }
464
465 Tc_(T), Tc_(bT); // XXX : ATTENTION : need to put back T in C
466
467 return 0; // FIXME : on suppose que tout OK
468#else
469 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
470 throw;
471#endif
472}
473
474int CoolProp_to_TRUST_generique::tppi_get_all_properties_h_IF97__(const MSpanD input, MLoiSpanD_h prop) const
475{
476#ifdef HAS_COOLPROP
477 const SpanD P = input.at("pression"), H = input.at("enthalpie");
478 const int sz = (int ) P.size();
479
480 bool has_prop = false;
481 for (auto &itr : prop)
482 if (itr.first == Loi_en_h::RHO || itr.first == Loi_en_h::T || itr.first == Loi_en_h::CP||
483 itr.first == Loi_en_h::MU || itr.first == Loi_en_h::LAMBDA || itr.first == Loi_en_h::BETA || itr.first == Loi_en_h::SIGMA)
484 {
485 has_prop = true;
486 break;
487 }
488
489 if (has_prop)
490 for (int i = 0; i < sz; i++)
491 {
492 fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i]); // SI units
493 for (auto &itr : prop)
494 {
495 Loi_en_h prop_ = itr.first;
496 SpanD span_ = itr.second;
497
498 if (prop_ == Loi_en_h::RHO) span_[i] = fluide->rhomass();
499 if (prop_ == Loi_en_h::T) span_[i] = fluide->T();
500 if (prop_ == Loi_en_h::CP) span_[i] = fluide->cpmass();
501 if (prop_ == Loi_en_h::MU) span_[i] = fluide->viscosity();
502 if (prop_ == Loi_en_h::LAMBDA) span_[i] = fluide->conductivity();
503 if (prop_ == Loi_en_h::BETA) span_[i] = fluide->isobaric_expansion_coefficient();
504 if (prop_ == Loi_en_h::SIGMA) span_[i] = fluide->surface_tension();
505 }
506 }
507
508 // derivees qui manquent ...
509 bool has_DP = false, has_DH = false;
510 for (auto &itr : prop)
511 if (itr.first == Loi_en_h::RHO_DP || itr.first == Loi_en_h::T_DP || itr.first == Loi_en_h::CP_DP ||
512 itr.first == Loi_en_h::MU_DP || itr.first == Loi_en_h::LAMBDA_DP || itr.first == Loi_en_h::SIGMA_DP)
513 {
514 has_DP = true;
515 break;
516 }
517
518 for (auto &itr : prop)
519 if (itr.first == Loi_en_h::RHO_DH || itr.first == Loi_en_h::T_DH || itr.first == Loi_en_h::CP_DH ||
520 itr.first == Loi_en_h::MU_DH || itr.first == Loi_en_h::LAMBDA_DH || itr.first == Loi_en_h::SIGMA_DH)
521 {
522 has_DH = true;
523 break;
524 }
525
526 if (has_DP)
527 for (int i = 0; i < sz; i++)
528 {
529 double plus_rho_ = EPS, plus_T_ = EPS, plus_cp_ = EPS, plus_mu_ = EPS, plus_lambda_ = EPS, plus_sigma_ = EPS;
530 fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i] * (1. + EPS));
531 for (auto &itr : prop)
532 {
533 Loi_en_h prop_ = itr.first;
534 if (prop_ == Loi_en_h::RHO_DP) plus_rho_ = fluide->rhomass();
535 if (prop_ == Loi_en_h::T_DP) plus_T_ = fluide->T();
536 if (prop_ == Loi_en_h::CP_DP) plus_cp_ = fluide->cpmass();
537 if (prop_ == Loi_en_h::MU_DP) plus_mu_ = fluide->viscosity();
538 if (prop_ == Loi_en_h::LAMBDA_DP) plus_lambda_ = fluide->conductivity();
539 if (prop_ == Loi_en_h::SIGMA_DP) plus_sigma_ = fluide->surface_tension();
540 }
541
542 double minus_rho_ = EPS, minus_T_ = EPS, minus_cp_ = EPS, minus_mu_ = EPS, minus_lambda_ = EPS, minus_sigma_ = EPS;
543 fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i] * (1. - EPS));
544 for (auto &itr : prop)
545 {
546 Loi_en_h prop_ = itr.first;
547 if (prop_ == Loi_en_h::RHO_DP) minus_rho_ = fluide->rhomass();
548 if (prop_ == Loi_en_h::T_DP) minus_T_ = fluide->T();
549 if (prop_ == Loi_en_h::CP_DP) minus_cp_ = fluide->cpmass();
550 if (prop_ == Loi_en_h::MU_DP) minus_mu_ = fluide->viscosity();
551 if (prop_ == Loi_en_h::LAMBDA_DP) minus_lambda_ = fluide->conductivity();
552 if (prop_ == Loi_en_h::SIGMA_DP) minus_sigma_ = fluide->surface_tension();
553 }
554
555 for (auto &itr : prop)
556 {
557 Loi_en_h prop_ = itr.first;
558 SpanD span_ = itr.second;
559
560 if (prop_ == Loi_en_h::RHO_DP) span_[i] = (plus_rho_ - minus_rho_) / ( 2 * EPS * P[i]);
561 if (prop_ == Loi_en_h::T_DP) span_[i] = (plus_T_ - minus_T_) / ( 2 * EPS * P[i]);
562 if (prop_ == Loi_en_h::CP_DP) span_[i] = (plus_cp_ - minus_cp_) / ( 2 * EPS * P[i]);
563 if (prop_ == Loi_en_h::MU_DP) span_[i] = (plus_mu_ - minus_mu_) / ( 2 * EPS * P[i]);
564 if (prop_ == Loi_en_h::LAMBDA_DP) span_[i] = (plus_lambda_ - minus_lambda_) / ( 2 * EPS * P[i]);
565 if (prop_ == Loi_en_h::SIGMA_DP) span_[i] = (plus_sigma_ - minus_sigma_) / ( 2 * EPS * P[i]);
566 }
567 }
568
569 if (has_DH)
570 for (int i = 0; i < sz; i++)
571 {
572 double plus_rho_ = EPS, plus_T_ = EPS, plus_cp_ = EPS, plus_mu_ = EPS, plus_lambda_ = EPS, plus_sigma_ = EPS;
573 fluide->update(CoolProp::HmassP_INPUTS, H[i] * (1. + EPS), P[i]);
574 for (auto &itr : prop)
575 {
576 Loi_en_h prop_ = itr.first;
577 if (prop_ == Loi_en_h::RHO_DH) plus_rho_ = fluide->rhomass();
578 if (prop_ == Loi_en_h::T_DH) plus_T_ = fluide->T();
579 if (prop_ == Loi_en_h::CP_DH) plus_cp_ = fluide->cpmass();
580 if (prop_ == Loi_en_h::MU_DH) plus_mu_ = fluide->viscosity();
581 if (prop_ == Loi_en_h::LAMBDA_DH) plus_lambda_ = fluide->conductivity();
582 if (prop_ == Loi_en_h::SIGMA_DH) plus_sigma_ = fluide->surface_tension();
583 }
584
585 double minus_rho_ = EPS, minus_T_ = EPS, minus_cp_ = EPS, minus_mu_ = EPS, minus_lambda_ = EPS, minus_sigma_ = EPS;
586 fluide->update(CoolProp::HmassP_INPUTS, H[i] * (1. - EPS), P[i]);
587 for (auto &itr : prop)
588 {
589 Loi_en_h prop_ = itr.first;
590 if (prop_ == Loi_en_h::RHO_DH) minus_rho_ = fluide->rhomass();
591 if (prop_ == Loi_en_h::T_DH) minus_T_ = fluide->T();
592 if (prop_ == Loi_en_h::CP_DH) minus_cp_ = fluide->cpmass();
593 if (prop_ == Loi_en_h::MU_DH) minus_mu_ = fluide->viscosity();
594 if (prop_ == Loi_en_h::LAMBDA_DH) minus_lambda_ = fluide->conductivity();
595 if (prop_ == Loi_en_h::SIGMA_DH) minus_sigma_ = fluide->surface_tension();
596 }
597
598 for (auto &itr : prop)
599 {
600 Loi_en_h prop_ = itr.first;
601 SpanD span_ = itr.second;
602
603 if (prop_ == Loi_en_h::RHO_DH) span_[i] = (plus_rho_ - minus_rho_) / ( 2 * EPS * H[i]);
604 if (prop_ == Loi_en_h::T_DH) span_[i] = (plus_T_ - minus_T_) / ( 2 * EPS * H[i]);
605 if (prop_ == Loi_en_h::CP_DH) span_[i] = (plus_cp_ - minus_cp_) / ( 2 * EPS * H[i]);
606 if (prop_ == Loi_en_h::MU_DH) span_[i] = (plus_mu_ - minus_mu_) / ( 2 * EPS * H[i]);
607 if (prop_ == Loi_en_h::LAMBDA_DH) span_[i] = (plus_lambda_ - minus_lambda_) / ( 2 * EPS * H[i]);
608 if (prop_ == Loi_en_h::SIGMA_DH) span_[i] = (plus_sigma_ - minus_sigma_) / ( 2 * EPS * H[i]);
609 }
610 }
611
612 return 0; // FIXME : on suppose que tout OK
613#else
614 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
615 throw;
616#endif
617}
618
619int CoolProp_to_TRUST_generique::tppi_get_all_properties_h__(const MSpanD input, MLoiSpanD_h prop) const
620{
621#ifdef HAS_COOLPROP
622 const SpanD P = input.at("pression"), H = input.at("enthalpie");
623 const int sz = (int ) P.size();
624
625 bool has_prop = false;
626 for (auto &itr : prop)
627 if (itr.first == Loi_en_h::RHO || itr.first == Loi_en_h::RHO_DP || itr.first == Loi_en_h::RHO_DH ||
628 itr.first == Loi_en_h::T || itr.first == Loi_en_h::T_DP || itr.first == Loi_en_h::T_DH ||
629 itr.first == Loi_en_h::CP || itr.first == Loi_en_h::CP_DP || itr.first == Loi_en_h::CP_DH ||
630 itr.first == Loi_en_h::MU || itr.first == Loi_en_h::LAMBDA || itr.first == Loi_en_h::BETA || itr.first == Loi_en_h::SIGMA)
631 {
632 has_prop = true;
633 break;
634 }
635
636 if (has_prop)
637 for (int i = 0; i < sz; i++)
638 {
639 fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i]); // SI units
640 for (auto &itr : prop)
641 {
642 Loi_en_h prop_ = itr.first;
643 SpanD span_ = itr.second;
644
645 if (prop_ == Loi_en_h::RHO) span_[i] = fluide->rhomass();
646 if (prop_ == Loi_en_h::RHO_DP) span_[i] = fluide->first_partial_deriv(CoolProp::iDmass, CoolProp::iP, CoolProp::iHmass);
647 if (prop_ == Loi_en_h::RHO_DH) span_[i] = fluide->first_partial_deriv(CoolProp::iDmass, CoolProp::iHmass, CoolProp::iP);
648 if (prop_ == Loi_en_h::T) span_[i] = fluide->T();
649 if (prop_ == Loi_en_h::T_DP) span_[i] = fluide->first_partial_deriv(CoolProp::iT, CoolProp::iP, CoolProp::iHmass);
650 if (prop_ == Loi_en_h::T_DH) span_[i] = fluide->first_partial_deriv(CoolProp::iT, CoolProp::iHmass, CoolProp::iP);
651 if (prop_ == Loi_en_h::CP) span_[i] = fluide->cpmass();
652 if (prop_ == Loi_en_h::CP_DP) span_[i] = fluide->first_partial_deriv(CoolProp::iCpmass, CoolProp::iP, CoolProp::iHmass);
653 if (prop_ == Loi_en_h::CP_DH) span_[i] = fluide->first_partial_deriv(CoolProp::iCpmass, CoolProp::iHmass, CoolProp::iP);
654 if (prop_ == Loi_en_h::MU) span_[i] = fluide->viscosity();
655 if (prop_ == Loi_en_h::LAMBDA) span_[i] = fluide->conductivity();
656 if (prop_ == Loi_en_h::BETA) span_[i] = fluide->isobaric_expansion_coefficient();
657 if (prop_ == Loi_en_h::SIGMA) span_[i] = fluide->surface_tension();
658 }
659 }
660
661 // derivees qui manquent ...
662 bool has_DP = false, has_DH = false;
663 for (auto &itr : prop)
664 if (itr.first == Loi_en_h::MU_DP || itr.first == Loi_en_h::LAMBDA_DP || itr.first == Loi_en_h::SIGMA_DP)
665 {
666 has_DP = true;
667 break;
668 }
669
670 for (auto &itr : prop)
671 if (itr.first == Loi_en_h::MU_DH || itr.first == Loi_en_h::LAMBDA_DH || itr.first == Loi_en_h::SIGMA_DH)
672 {
673 has_DH = true;
674 break;
675 }
676
677 if (has_DP)
678 for (int i = 0; i < sz; i++)
679 {
680 double plus_mu_ = EPS, plus_lambda_ = EPS, plus_sigma_ = EPS;
681 fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i] * (1. + EPS));
682 for (auto &itr : prop)
683 {
684 Loi_en_h prop_ = itr.first;
685 if (prop_ == Loi_en_h::MU_DP) plus_mu_ = fluide->viscosity();
686 if (prop_ == Loi_en_h::LAMBDA_DP) plus_lambda_ = fluide->conductivity();
687 if (prop_ == Loi_en_h::SIGMA_DP) plus_sigma_ = fluide->surface_tension();
688 }
689
690 double minus_mu_ = EPS, minus_lambda_ = EPS, minus_sigma_ = EPS;
691 fluide->update(CoolProp::HmassP_INPUTS, H[i], P[i] * (1. - EPS));
692 for (auto &itr : prop)
693 {
694 Loi_en_h prop_ = itr.first;
695 if (prop_ == Loi_en_h::MU_DP) minus_mu_ = fluide->viscosity();
696 if (prop_ == Loi_en_h::LAMBDA_DP) minus_lambda_ = fluide->conductivity();
697 if (prop_ == Loi_en_h::SIGMA_DP) minus_sigma_ = fluide->surface_tension();
698 }
699
700 for (auto &itr : prop)
701 {
702 Loi_en_h prop_ = itr.first;
703 SpanD span_ = itr.second;
704
705 if (prop_ == Loi_en_h::MU_DP) span_[i] = (plus_mu_ - minus_mu_) / ( 2 * EPS * P[i]);
706 if (prop_ == Loi_en_h::LAMBDA_DP) span_[i] = (plus_lambda_ - minus_lambda_) / ( 2 * EPS * P[i]);
707 if (prop_ == Loi_en_h::SIGMA_DP) span_[i] = (plus_sigma_ - minus_sigma_) / ( 2 * EPS * P[i]);
708 }
709 }
710
711 if (has_DH)
712 for (int i = 0; i < sz; i++)
713 {
714 double plus_mu_ = EPS, plus_lambda_ = EPS, plus_sigma_ = EPS;
715 fluide->update(CoolProp::HmassP_INPUTS, H[i] * (1. + EPS), P[i]);
716 for (auto &itr : prop)
717 {
718 Loi_en_h prop_ = itr.first;
719 if (prop_ == Loi_en_h::MU_DH) plus_mu_ = fluide->viscosity();
720 if (prop_ == Loi_en_h::LAMBDA_DH) plus_lambda_ = fluide->conductivity();
721 if (prop_ == Loi_en_h::SIGMA_DH) plus_sigma_ = fluide->surface_tension();
722 }
723
724 double minus_mu_ = EPS, minus_lambda_ = EPS, minus_sigma_ = EPS;
725 fluide->update(CoolProp::HmassP_INPUTS, H[i] * (1. - EPS), P[i]);
726 for (auto &itr : prop)
727 {
728 Loi_en_h prop_ = itr.first;
729 if (prop_ == Loi_en_h::MU_DH) minus_mu_ = fluide->viscosity();
730 if (prop_ == Loi_en_h::LAMBDA_DH) minus_lambda_ = fluide->conductivity();
731 if (prop_ == Loi_en_h::SIGMA_DH) minus_sigma_ = fluide->surface_tension();
732 }
733
734 for (auto &itr : prop)
735 {
736 Loi_en_h prop_ = itr.first;
737 SpanD span_ = itr.second;
738
739 if (prop_ == Loi_en_h::MU_DH) span_[i] = (plus_mu_ - minus_mu_) / ( 2 * EPS * H[i]);
740 if (prop_ == Loi_en_h::LAMBDA_DH) span_[i] = (plus_lambda_ - minus_lambda_) / ( 2 * EPS * H[i]);
741 if (prop_ == Loi_en_h::SIGMA_DH) span_[i] = (plus_sigma_ - minus_sigma_) / ( 2 * EPS * H[i]);
742 }
743 }
744
745 return 0; // FIXME : on suppose que tout OK
746#else
747 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
748 throw;
749#endif
750}
751
752// methodes particulieres par application pour gagner en performance : utilisees dans F5 (pour le moment !)
753int CoolProp_to_TRUST_generique::tppi_get_all_prop_loi_F5(const MSpanD input, MLoiSpanD_h inter, int ncomp, int ind, bool is_liq) const
754{
755#ifdef HAS_COOLPROP
756 const SpanD P = input.at("pression"), H = is_liq ? input.at("H_L") : input.at("H_V");
757
758#ifndef NDEBUG
759 assert(ncomp == 1);
760 for (auto &itr : inter) assert((int )H.size() == ncomp * (int )itr.second.size());
761#endif
762
763 MSpanD input_ = { {"pression", P }, { "enthalpie", H} };
764
765 if (fluide->backend_name() == "IF97Backend")
766 return tppi_get_all_properties_h_IF97__(input_, inter);
767 else
768 return tppi_get_all_properties_h__(input_, inter);
769#else
770 Cerr << "CoolProp_to_TRUST_generique::" << __func__ << " should not be called since TRUST is not compiled with the CoolProp library !!! " << finl;
771 throw;
772#endif
773}
CoolProp_to_TRUST_generique::set_fluide_generique
void set_fluide_generique(const char *const model_name, const char *const fluid_name) override
Definition CoolProp_to_TRUST_generique.cpp:24
CoolProp_to_TRUST_generique::tppi_get_all_prop_loi_F5
int tppi_get_all_prop_loi_F5(const MSpanD, MLoiSpanD_h, int ncomp=1, int id=0, bool is_liq=true) const override
Definition CoolProp_to_TRUST_generique.cpp:753
CoolProp_to_TRUST_generique::tppi_get_all_pb_multiphase_pT
int tppi_get_all_pb_multiphase_pT(const MSpanD, MLoiSpanD, MLoiSpanD, int ncomp=1, int id=0) const override
Definition CoolProp_to_TRUST_generique.cpp:415
CoolProp_to_TRUST_generique::tppi_get_CPMLB_pb_multiphase_pT
int tppi_get_CPMLB_pb_multiphase_pT(const MSpanD, MLoiSpanD, int ncomp=1, int id=0) const override
Definition CoolProp_to_TRUST_generique.cpp:381
Process::exit
static void exit(int exit_code=-1)
Routine de sortie de TRUST dans une region Kokkos.
Definition Process.cpp:455