CoolProp 6.8.1dev
An open-source fluid property and humid air property database
CubicBackend.h
Go to the documentation of this file.
1/*
2
3The goal of this backend is to allow the Helmholtz-based structure for cubics and to replace the entire
4multi-fluid model with a one-fluid model. The other changes are relatively trivial. The primary
5change is to replace the core residual Helmholtz energy derivatives from HelmholtzEOSMixtureBackend
6with the derivatives from this class.
7
8The core code for the Helmholtz translations is from the publication
9"Helmholtz energy translations for common cubic equations of state for use in one-fluid and multi-fluid mixture models"
10by Ian H. Bell and Andreas Jaeger, J. Res. NIST, 2016
11
12*/
13
14#ifndef CUBICBACKEND_H_
15#define CUBICBACKEND_H_
16
17#include "CoolPropTools.h"
18#include "DataStructures.h"
19#include "GeneralizedCubic.h"
20#include "CubicsLibrary.h"
21#include "Configuration.h"
22#include "AbstractState.h"
23#include "Backends/Helmholtz/HelmholtzEOSMixtureBackend.h"
24#include "Exceptions.h"
25#include <vector>
26
27namespace CoolProp {
28
29// Forward declaration for use in initialization of AbstractCubicBackend
30class CubicResidualHelmholtz;
31
32class AbstractCubicBackend : public HelmholtzEOSMixtureBackend
33{
34 protected:
35 shared_ptr<AbstractCubic> cubic;
36 std::vector<CubicLibrary::CubicsValues> components;
37 public:
39 void setup(bool generate_SatL_and_SatV = true);
40
42 void set_alpha_from_components();
43
45 void set_alpha0_from_components();
46
48 shared_ptr<AbstractCubic>& get_cubic() {
49 return cubic;
50 };
51
52 std::vector<std::string> calc_fluid_names(void);
53
54 bool using_mole_fractions(void) {
55 return true;
56 };
57 bool using_mass_fractions(void) {
58 return false;
59 };
60 bool using_volu_fractions(void) {
61 return false;
62 };
63
64 void set_mass_fractions(const std::vector<CoolPropDbl>& mass_fractions) {
65 throw NotImplementedError("Mass composition has not been implemented.");
66 };
67 void set_volu_fractions(const std::vector<CoolPropDbl>& volu_fractions) {
68 throw NotImplementedError("Volume composition has not been implemented.");
69 };
70 const std::vector<CoolPropDbl>& get_mole_fractions(void) {
71 return this->mole_fractions;
72 };
73
74 const double get_fluid_constant(std::size_t i, parameters param) const {
75 switch (param) {
76 case iP_critical:
77 return cubic->get_pc()[i];
78 case iT_reducing:
79 case iT_critical:
80 return cubic->get_Tc()[i];
81 case iacentric_factor:
82 return cubic->get_acentric()[i];
83 case imolar_mass:
84 return components[i].molemass;
85 case iT_triple:
86 return HelmholtzEOSMixtureBackend::get_components()[i].EOS().sat_min_liquid.T; // From the base class data structure
87 case iP_triple:
88 return HelmholtzEOSMixtureBackend::get_components()[i].EOS().sat_min_liquid.p; // From the base class data structure
89 case irhomolar_reducing:
90 case irhomolar_critical:
91 return components[i].rhomolarc;
92 case igas_constant:
93 return get_config_double(R_U_CODATA);
94 default:
95 throw ValueError(format("I don't know what to do with this fluid constant: %s", get_parameter_information(param, "short").c_str()));
96 }
97 }
98
100 std::string fluid_param_string(const std::string&);
102 CoolPropDbl calc_gas_constant(void) {
103 return cubic->get_R_u();
104 };
106 SimpleState calc_reducing_state_nocache(const std::vector<CoolPropDbl>& mole_fractions) {
107 SimpleState reducing;
108 reducing.T = cubic->get_Tr();
109 reducing.rhomolar = cubic->get_rhor();
110 return reducing;
111 };
112 CoolPropDbl calc_reduced_density(void) {
113 return _rhomolar / get_cubic()->get_rhor();
114 };
115 CoolPropDbl calc_reciprocal_reduced_temperature(void) {
116 return get_cubic()->get_Tr() / _T;
117 };
118 std::vector<double> spinodal_densities();
119
120 CoolPropDbl calc_T_critical(void) {
121 if (is_pure_or_pseudopure) {
122 return cubic->get_Tc()[0];
123 } else {
124 return HelmholtzEOSMixtureBackend::calc_T_critical();
125 }
126 };
127 CoolPropDbl calc_p_critical(void) {
128 if (is_pure_or_pseudopure) {
129 return cubic->get_pc()[0];
130 } else {
131 return HelmholtzEOSMixtureBackend::calc_p_critical();
132 }
133 };
134 CoolPropDbl calc_acentric_factor(void) {
135 if (is_pure_or_pseudopure) {
136 return cubic->get_acentric()[0];
137 } else {
138 throw ValueError("acentric factor cannot be calculated for mixtures");
139 }
140 }
141 CoolPropDbl calc_rhomolar_critical(void) {
142 if (is_pure_or_pseudopure) {
143 // Curve fit from all the pure fluids in CoolProp (thanks to recommendation of A. Kazakov)
144 double v_c_Lmol = 2.14107171795 * (cubic->get_Tc()[0] / cubic->get_pc()[0] * 1000) + 0.00773144012514; // [L/mol]
145 return 1 / (v_c_Lmol / 1000.0);
146 } else {
147 return HelmholtzEOSMixtureBackend::calc_rhomolar_critical();
148 }
149 };
150
153 void get_linear_reducing_parameters(double& rhomolar, double& T);
154
156 void get_critical_point_starting_values(double& delta0, double& tau0);
157
159 void get_critical_point_search_radii(double& R_delta, double& R_tau);
160
162 bool get_critical_is_terminated(double& delta, double& tau);
163
164 CoolPropDbl calc_alphar_deriv_nocache(const int nTau, const int nDelta, const std::vector<CoolPropDbl>& mole_fractions, const CoolPropDbl& tau,
165 const CoolPropDbl& delta);
166
168 CoolPropDbl calc_pressure_nocache(CoolPropDbl T, CoolPropDbl rhomolar);
169
171 virtual void update_DmolarT();
172
173 virtual void update(CoolProp::input_pairs input_pair, double value1, double value2);
174
197 void rho_Tp_cubic(CoolPropDbl T, CoolPropDbl p, int& Nsolns, double& rho0, double& rho1, double& rho2);
198
200 CoolPropDbl solver_rho_Tp_SRK(CoolPropDbl T, CoolPropDbl p, phases phase) {
201 return solver_rho_Tp(T, p);
202 };
208 CoolPropDbl solver_rho_Tp(CoolPropDbl T, CoolPropDbl p, CoolPropDbl rho_guess = -1);
209
210 CoolPropDbl solver_rho_Tp_global(CoolPropDbl T, CoolPropDbl p, CoolPropDbl rhomax);
211
213 void update_TPD_state() {
214 TPD_state.reset(get_copy());
215 };
216
218 void saturation(CoolProp::input_pairs inputs);
219
220 CoolPropDbl calc_molar_mass(void);
221
222 void set_binary_interaction_double(const std::size_t i1, const std::size_t i2, const std::string& parameter, const double value);
223 double get_binary_interaction_double(const std::size_t i1, const std::size_t i2, const std::string& parameter);
224
225 void set_binary_interaction_double(const std::string& CAS1, const std::string& CAS2, const std::string& parameter, const double value) {
226 throw ValueError("set_binary_interaction_double not defined for AbstractCubic not defined for CAS #");
227 }
228 double get_binary_interaction_double(const std::string& CAS1, const std::string& CAS2, const std::string& parameter) {
229 throw ValueError("get_binary_interaction_double not defined for AbstractCubic not defined for CAS #");
230 };
231
232 // Return a 1-1 copy of this class
233 virtual HelmholtzEOSMixtureBackend* get_copy(bool generate_SatL_and_SatV = true) = 0;
234
235 // Copy the entire kij matrix from another instance in one shot
236 void copy_k(AbstractCubicBackend* donor);
237
238 //
239 void copy_all_alpha_functions(AbstractCubicBackend* donor);
240
242 void copy_internals(AbstractCubicBackend& donor);
243
244 // Set the cubic alpha function's constants:
245 void set_cubic_alpha_C(const size_t i, const std::string& parameter, const double c1, const double c2, const double c3);
246
247 // Set fluid parameter (currently the volume translation parameter)
248 void set_fluid_parameter_double(const size_t i, const std::string& parameter, const double value);
249
250 // Get fluid parameter (currently the volume translation parameter)
251 double get_fluid_parameter_double(const size_t i, const std::string& parameter);
252};
253
254class SRKBackend : public AbstractCubicBackend
255{
256
257 public:
258 SRKBackend(const std::vector<double>& Tc, const std::vector<double>& pc, const std::vector<double>& acentric, double R_u,
259 bool generate_SatL_and_SatV = true) {
260 cubic.reset(new SRK(Tc, pc, acentric, R_u));
261 setup(generate_SatL_and_SatV);
262 };
263 SRKBackend(double Tc, double pc, double acentric, double R_u, bool generate_SatL_and_SatV = true) {
264 cubic.reset(new SRK(Tc, pc, acentric, R_u));
265 setup(generate_SatL_and_SatV);
266 }
267 SRKBackend(const std::vector<std::string> fluid_identifiers, const double R_u = get_config_double(R_U_CODATA),
268 bool generate_SatL_and_SatV = true) {
269 std::vector<double> Tc, pc, acentric;
270 N = fluid_identifiers.size();
271 components.resize(N);
272 for (std::size_t i = 0; i < fluid_identifiers.size(); ++i) {
273 components[i] = CubicLibrary::get_cubic_values(fluid_identifiers[i]);
274 Tc.push_back(components[i].Tc);
275 pc.push_back(components[i].pc);
276 acentric.push_back(components[i].acentric);
277 }
278 cubic.reset(new SRK(Tc, pc, acentric, R_u));
279 setup(generate_SatL_and_SatV);
280 }
281 HelmholtzEOSMixtureBackend* get_copy(bool generate_SatL_and_SatV = true) {
282 AbstractCubicBackend* ACB = new SRKBackend(cubic->get_Tc(), cubic->get_pc(), cubic->get_acentric(), cubic->get_R_u(), generate_SatL_and_SatV);
283 ACB->copy_internals(*this);
284 return static_cast<HelmholtzEOSMixtureBackend*>(ACB);
285 }
286 std::string backend_name(void) {
287 return get_backend_string(SRK_BACKEND);
288 }
289};
290
291class PengRobinsonBackend : public AbstractCubicBackend
292{
293
294 public:
295 PengRobinsonBackend(){}; // Default constructor (make sure you know what you are doing)
296 PengRobinsonBackend(const std::vector<double>& Tc, const std::vector<double>& pc, const std::vector<double>& acentric, double R_u,
297 bool generate_SatL_and_SatV = true) {
298 cubic.reset(new PengRobinson(Tc, pc, acentric, R_u));
299 setup(generate_SatL_and_SatV);
300 };
301 PengRobinsonBackend(double Tc, double pc, double acentric, double R_u, bool generate_SatL_and_SatV = true) {
302 cubic.reset(new PengRobinson(Tc, pc, acentric, R_u));
303 setup(generate_SatL_and_SatV);
304 };
305 PengRobinsonBackend(const std::vector<std::string> fluid_identifiers, const double R_u = get_config_double(R_U_CODATA),
306 bool generate_SatL_and_SatV = true) {
307 std::vector<double> Tc, pc, acentric;
308 N = fluid_identifiers.size();
309 components.resize(N);
310 for (std::size_t i = 0; i < fluid_identifiers.size(); ++i) {
311 components[i] = CubicLibrary::get_cubic_values(fluid_identifiers[i]);
312 Tc.push_back(components[i].Tc);
313 pc.push_back(components[i].pc);
314 acentric.push_back(components[i].acentric);
315 }
316 cubic.reset(new PengRobinson(Tc, pc, acentric, R_u));
317 setup(generate_SatL_and_SatV);
318 };
319 HelmholtzEOSMixtureBackend* get_copy(bool generate_SatL_and_SatV = true) {
320 AbstractCubicBackend* ACB =
321 new PengRobinsonBackend(cubic->get_Tc(), cubic->get_pc(), cubic->get_acentric(), cubic->get_R_u(), generate_SatL_and_SatV);
322 ACB->copy_internals(*this);
323 return static_cast<HelmholtzEOSMixtureBackend*>(ACB);
324 }
325 std::string backend_name(void) {
326 return get_backend_string(PR_BACKEND);
327 }
328};
329
334class CubicResidualHelmholtz : public ResidualHelmholtz
335{
336
337 protected:
338 AbstractCubicBackend* ACB;
339
340 public:
341 CubicResidualHelmholtz() {
342 ACB = NULL;
343 };
344 CubicResidualHelmholtz(AbstractCubicBackend* ACB) : ACB(ACB){};
345
346 // copy assignment
347 CubicResidualHelmholtz& operator=(CubicResidualHelmholtz& other) {
348 ACB = other.ACB;
349 return *this;
350 }
351
353 virtual HelmholtzDerivatives all(HelmholtzEOSMixtureBackend& HEOS, const std::vector<CoolPropDbl>& mole_fractions, double tau, double delta,
354 bool cache_values = false) {
355 HelmholtzDerivatives a;
356 std::vector<double> z = std::vector<double>(mole_fractions.begin(), mole_fractions.end());
357 shared_ptr<AbstractCubic>& cubic = ACB->get_cubic();
358 a.alphar = cubic->alphar(tau, delta, z, 0, 0);
359 a.dalphar_dtau = cubic->alphar(tau, delta, z, 1, 0);
360 a.dalphar_ddelta = cubic->alphar(tau, delta, z, 0, 1);
361 a.d2alphar_dtau2 = cubic->alphar(tau, delta, z, 2, 0);
362 a.d2alphar_ddelta_dtau = cubic->alphar(tau, delta, z, 1, 1);
363 a.d2alphar_ddelta2 = cubic->alphar(tau, delta, z, 0, 2);
364 a.d3alphar_dtau3 = cubic->alphar(tau, delta, z, 3, 0);
365 a.d3alphar_ddelta_dtau2 = cubic->alphar(tau, delta, z, 2, 1);
366 a.d3alphar_ddelta2_dtau = cubic->alphar(tau, delta, z, 1, 2);
367 a.d3alphar_ddelta3 = cubic->alphar(tau, delta, z, 0, 3);
368 a.d4alphar_dtau4 = cubic->alphar(tau, delta, z, 4, 0);
369 a.d4alphar_ddelta_dtau3 = cubic->alphar(tau, delta, z, 3, 1);
370 a.d4alphar_ddelta2_dtau2 = cubic->alphar(tau, delta, z, 2, 2);
371 a.d4alphar_ddelta3_dtau = cubic->alphar(tau, delta, z, 1, 3);
372 a.d4alphar_ddelta4 = cubic->alphar(tau, delta, z, 0, 4);
373 return a;
374 }
375 virtual CoolPropDbl dalphar_dxi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
376 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 0, i, xN_flag == XN_INDEPENDENT);
377 }
378 virtual CoolPropDbl d2alphar_dxi_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
379 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 1, 0, i, xN_flag == XN_INDEPENDENT);
380 }
381 virtual CoolPropDbl d2alphar_dxi_dDelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
382 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 1, i, xN_flag == XN_INDEPENDENT);
383 }
384 virtual CoolPropDbl d3alphar_dxi_dTau2(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
385 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 2, 0, i, xN_flag == XN_INDEPENDENT);
386 }
387 virtual CoolPropDbl d3alphar_dxi_dDelta_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
388 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 1, 1, i, xN_flag == XN_INDEPENDENT);
389 }
390 virtual CoolPropDbl d3alphar_dxi_dDelta2(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
391 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 2, i, xN_flag == XN_INDEPENDENT);
392 }
393
394 virtual CoolPropDbl d2alphardxidxj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
395 return ACB->get_cubic()->d2_alphar_dxidxj(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 0, i, j,
396 xN_flag == XN_INDEPENDENT);
397 }
398 virtual CoolPropDbl d3alphar_dxi_dxj_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
399 return ACB->get_cubic()->d2_alphar_dxidxj(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 1, 0, i, j,
400 xN_flag == XN_INDEPENDENT);
401 }
402 virtual CoolPropDbl d3alphar_dxi_dxj_dDelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
403 return ACB->get_cubic()->d2_alphar_dxidxj(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 1, i, j,
404 xN_flag == XN_INDEPENDENT);
405 }
406
407 virtual CoolPropDbl d4alphar_dxi_dTau3(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
408 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 3, 0, i, xN_flag == XN_INDEPENDENT);
409 }
410 virtual CoolPropDbl d4alphar_dxi_dDelta2_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
411 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 1, 2, i, xN_flag == XN_INDEPENDENT);
412 }
413 virtual CoolPropDbl d4alphar_dxi_dDelta_dTau2(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
414 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 2, 1, i, xN_flag == XN_INDEPENDENT);
415 }
416 virtual CoolPropDbl d4alphar_dxi_dDelta3(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
417 return ACB->get_cubic()->d_alphar_dxi(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 3, i, xN_flag == XN_INDEPENDENT);
418 }
419 virtual CoolPropDbl d4alphar_dxi_dxj_dTau2(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
420 return ACB->get_cubic()->d2_alphar_dxidxj(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 2, 0, i, j,
421 xN_flag == XN_INDEPENDENT);
422 }
423 virtual CoolPropDbl d4alphar_dxi_dxj_dDelta_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
424 return ACB->get_cubic()->d2_alphar_dxidxj(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 1, 1, i, j,
425 xN_flag == XN_INDEPENDENT);
426 }
427 virtual CoolPropDbl d4alphar_dxi_dxj_dDelta2(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
428 return ACB->get_cubic()->d2_alphar_dxidxj(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 2, i, j,
429 xN_flag == XN_INDEPENDENT);
430 }
431 virtual CoolPropDbl d3alphardxidxjdxk(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, std::size_t k,
432 x_N_dependency_flag xN_flag) {
433 return ACB->get_cubic()->d3_alphar_dxidxjdxk(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 0, i, j, k,
434 xN_flag == XN_INDEPENDENT);
435 }
436 virtual CoolPropDbl d4alphar_dxi_dxj_dxk_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, std::size_t k,
437 x_N_dependency_flag xN_flag) {
438 return ACB->get_cubic()->d3_alphar_dxidxjdxk(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 1, 0, i, j, k,
439 xN_flag == XN_INDEPENDENT);
440 }
441 virtual CoolPropDbl d4alphar_dxi_dxj_dxk_dDelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, std::size_t k,
442 x_N_dependency_flag xN_flag) {
443 return ACB->get_cubic()->d3_alphar_dxidxjdxk(HEOS.tau(), HEOS.delta(), HEOS.get_mole_fractions_doubleref(), 0, 1, i, j, k,
444 xN_flag == XN_INDEPENDENT);
445 }
446};
447
448} /* namespace CoolProp */
449#endif /* CUBICBACKEND_H_ */
Generated on Wed Jun 4 2025 04:12:26 for CoolProp by doxygen 1.9.4