CoolProp  6.6.0
An open-source fluid property and humid air property database
CubicBackend.h
Go to the documentation of this file.
1 /*
2 
3 The goal of this backend is to allow the Helmholtz-based structure for cubics and to replace the entire
4 multi-fluid model with a one-fluid model. The other changes are relatively trivial. The primary
5 change is to replace the core residual Helmholtz energy derivatives from HelmholtzEOSMixtureBackend
6 with the derivatives from this class.
7 
8 The 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"
10 by 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 
27 namespace CoolProp {
28 
29 // Forward declaration for use in initialization of AbstractCubicBackend
30 class CubicResidualHelmholtz;
31 
32 class 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 
254 class 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 
291 class 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 
334 class 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 Nov 29 2023 14:26:23 for CoolProp by doxygen 1.9.1