CoolProp
AbstractState.h
1 /*
2  * AbstractState.h
3  *
4  * Created on: 21 Dec 2013
5  * Author: jowr
6  */
7 
8 #ifndef ABSTRACTSTATE_H_
9 #define ABSTRACTSTATE_H_
10 
11 #include "CachedElement.h"
12 #include "Exceptions.h"
13 #include "DataStructures.h"
14 #include "PhaseEnvelope.h"
15 #include "crossplatform_shared_ptr.h"
16 
17 #include <numeric>
18 
19 namespace CoolProp {
20 
24 public:
25  std::vector<double> tau,
26  delta,
27  M1;
28 };
29 
33 public:
34  double T,
35  p,
36  rhomolar,
37  hmolar,
38  smolar,
39  rhomolar_liq,
40  rhomolar_vap;
41  std::vector<double> x,
42  y;
44  clear();
45  };
46  void clear() {
47  T = _HUGE; p = _HUGE; rhomolar = _HUGE; hmolar = _HUGE; smolar = _HUGE;
48  rhomolar_liq = _HUGE; rhomolar_vap = _HUGE; x.clear(), y.clear();
49  }
50 };
51 
53 
71 protected:
72 
77 
78  bool isSupercriticalPhase(void){
79  return (this->_phase == iphase_supercritical || this->_phase == iphase_supercritical_liquid || this->_phase == iphase_supercritical_gas);
80  }
81 
82  bool isHomogeneousPhase(void){
83  return (this->_phase == iphase_liquid || this->_phase == iphase_gas || isSupercriticalPhase() || this->_phase == iphase_critical_point);
84  }
85 
86  bool isTwoPhase(void){
87  return (this->_phase == iphase_twophase);
88  }
89 
91  SimpleState _critical, _reducing;
92 
95 
98 
100  double _rhomolar, _T, _p, _Q, _R;
101 
102  CachedElement _tau, _delta;
103 
105  CachedElement _viscosity, _conductivity, _surface_tension;
106 
107  CachedElement _hmolar, _smolar, _umolar, _logp, _logrhomolar, _cpmolar, _cp0molar, _cvmolar, _speed_sound, _gibbsmolar, _helmholtzmolar;
108 
110  CachedElement _hmolar_residual, _smolar_residual, _gibbsmolar_residual;
111 
113  CachedElement _hmolar_excess, _smolar_excess, _gibbsmolar_excess, _umolar_excess, _volumemolar_excess, _helmholtzmolar_excess;
114 
116  CachedElement _rhoLanc, _rhoVanc, _pLanc, _pVanc, _TLanc, _TVanc;
117 
118  CachedElement _fugacity_coefficient;
119 
121  CachedElement _rho_spline, _drho_spline_dh__constp, _drho_spline_dp__consth;
122 
123 
125  CachedElement _alpha0, _dalpha0_dTau, _dalpha0_dDelta, _d2alpha0_dTau2, _d2alpha0_dDelta_dTau,
126  _d2alpha0_dDelta2, _d3alpha0_dTau3, _d3alpha0_dDelta_dTau2, _d3alpha0_dDelta2_dTau,
127  _d3alpha0_dDelta3, _alphar, _dalphar_dTau, _dalphar_dDelta, _d2alphar_dTau2, _d2alphar_dDelta_dTau,
128  _d2alphar_dDelta2, _d3alphar_dTau3, _d3alphar_dDelta_dTau2, _d3alphar_dDelta2_dTau,
129  _d3alphar_dDelta3, _d4alphar_dTau4, _d4alphar_dDelta_dTau3, _d4alphar_dDelta2_dTau2,
130  _d4alphar_dDelta3_dTau, _d4alphar_dDelta4;
131 
132  CachedElement _dalphar_dDelta_lim, _d2alphar_dDelta2_lim,
133  _d2alphar_dDelta_dTau_lim, _d3alphar_dDelta2_dTau_lim;
134 
137 
138  // ----------------------------------------
139  // Property accessors to be optionally implemented by the backend
140  // for properties that are not always calculated
141  // ----------------------------------------
143  virtual CoolPropDbl calc_hmolar(void){ throw NotImplementedError("calc_hmolar is not implemented for this backend"); };
145  virtual CoolPropDbl calc_hmolar_residual(void){ throw NotImplementedError("calc_hmolar_residual is not implemented for this backend"); };
147  virtual CoolPropDbl calc_smolar(void){ throw NotImplementedError("calc_smolar is not implemented for this backend"); };
149  virtual CoolPropDbl calc_smolar_residual(void){ throw NotImplementedError("calc_smolar_residual is not implemented for this backend"); };
151  virtual CoolPropDbl calc_umolar(void){ throw NotImplementedError("calc_umolar is not implemented for this backend"); };
153  virtual CoolPropDbl calc_cpmolar(void){ throw NotImplementedError("calc_cpmolar is not implemented for this backend"); };
155  virtual CoolPropDbl calc_cpmolar_idealgas(void){ throw NotImplementedError("calc_cpmolar_idealgas is not implemented for this backend"); };
157  virtual CoolPropDbl calc_cvmolar(void){ throw NotImplementedError("calc_cvmolar is not implemented for this backend"); };
159  virtual CoolPropDbl calc_gibbsmolar(void){ throw NotImplementedError("calc_gibbsmolar is not implemented for this backend"); };
161  virtual CoolPropDbl calc_gibbsmolar_residual(void){ throw NotImplementedError("calc_gibbsmolar_residual is not implemented for this backend"); };
163  virtual CoolPropDbl calc_helmholtzmolar(void){ throw NotImplementedError("calc_helmholtzmolar is not implemented for this backend"); };
165  virtual CoolPropDbl calc_speed_sound(void){ throw NotImplementedError("calc_speed_sound is not implemented for this backend"); };
167  virtual CoolPropDbl calc_isothermal_compressibility(void){ throw NotImplementedError("calc_isothermal_compressibility is not implemented for this backend"); };
169  virtual CoolPropDbl calc_isobaric_expansion_coefficient(void){ throw NotImplementedError("calc_isobaric_expansion_coefficient is not implemented for this backend"); };
171  virtual CoolPropDbl calc_isentropic_expansion_coefficient(void) { throw NotImplementedError("calc_isentropic_expansion_coefficient is not implemented for this backend"); };
173  virtual CoolPropDbl calc_viscosity(void){ throw NotImplementedError("calc_viscosity is not implemented for this backend"); };
175  virtual CoolPropDbl calc_conductivity(void){ throw NotImplementedError("calc_conductivity is not implemented for this backend"); };
177  virtual CoolPropDbl calc_surface_tension(void){ throw NotImplementedError("calc_surface_tension is not implemented for this backend"); };
179  virtual CoolPropDbl calc_molar_mass(void){ throw NotImplementedError("calc_molar_mass is not implemented for this backend"); };
181  virtual CoolPropDbl calc_acentric_factor(void){ throw NotImplementedError("calc_acentric_factor is not implemented for this backend"); };
183  virtual CoolPropDbl calc_pressure(void){ throw NotImplementedError("calc_pressure is not implemented for this backend"); };
185  virtual CoolPropDbl calc_gas_constant(void){ throw NotImplementedError("calc_gas_constant is not implemented for this backend"); };
187  virtual CoolPropDbl calc_fugacity_coefficient(std::size_t i){ throw NotImplementedError("calc_fugacity_coefficient is not implemented for this backend"); };
189  virtual std::vector<CoolPropDbl> calc_fugacity_coefficients(){ throw NotImplementedError("calc_fugacity_coefficients is not implemented for this backend"); };
191  virtual CoolPropDbl calc_fugacity(std::size_t i){ throw NotImplementedError("calc_fugacity is not implemented for this backend"); };
193  virtual CoolPropDbl calc_chemical_potential(std::size_t i) { throw NotImplementedError("calc_chemical_potential is not implemented for this backend"); };
195  virtual CoolPropDbl calc_PIP(void){ throw NotImplementedError("calc_PIP is not implemented for this backend"); };
196 
197  // Excess properties
199  virtual void calc_excess_properties(void) { throw NotImplementedError("calc_excess_properties is not implemented for this backend"); };
200 
201  // Derivatives of residual helmholtz energy
203  virtual CoolPropDbl calc_alphar(void){ throw NotImplementedError("calc_alphar is not implemented for this backend"); };
205  virtual CoolPropDbl calc_dalphar_dDelta(void){ throw NotImplementedError("calc_dalphar_dDelta is not implemented for this backend"); };
207  virtual CoolPropDbl calc_dalphar_dTau(void){ throw NotImplementedError("calc_dalphar_dTau is not implemented for this backend"); };
209  virtual CoolPropDbl calc_d2alphar_dDelta2(void){ throw NotImplementedError("calc_d2alphar_dDelta2 is not implemented for this backend"); };
211  virtual CoolPropDbl calc_d2alphar_dDelta_dTau(void){ throw NotImplementedError("calc_d2alphar_dDelta_dTau is not implemented for this backend"); };
213  virtual CoolPropDbl calc_d2alphar_dTau2(void){ throw NotImplementedError("calc_d2alphar_dTau2 is not implemented for this backend"); };
215  virtual CoolPropDbl calc_d3alphar_dDelta3(void){ throw NotImplementedError("calc_d3alphar_dDelta3 is not implemented for this backend"); };
217  virtual CoolPropDbl calc_d3alphar_dDelta2_dTau(void){ throw NotImplementedError("calc_d3alphar_dDelta2_dTau is not implemented for this backend"); };
219  virtual CoolPropDbl calc_d3alphar_dDelta_dTau2(void){ throw NotImplementedError("calc_d3alphar_dDelta_dTau2 is not implemented for this backend"); };
221  virtual CoolPropDbl calc_d3alphar_dTau3(void){ throw NotImplementedError("calc_d3alphar_dTau3 is not implemented for this backend"); };
222 
224  virtual CoolPropDbl calc_d4alphar_dDelta4(void){ throw NotImplementedError("calc_d4alphar_dDelta4 is not implemented for this backend"); };
226  virtual CoolPropDbl calc_d4alphar_dDelta3_dTau(void){ throw NotImplementedError("calc_d4alphar_dDelta3_dTau is not implemented for this backend"); };
228  virtual CoolPropDbl calc_d4alphar_dDelta2_dTau2(void){ throw NotImplementedError("calc_d4alphar_dDelta2_dTau2 is not implemented for this backend"); };
230  virtual CoolPropDbl calc_d4alphar_dDelta_dTau3(void){ throw NotImplementedError("calc_d4alphar_dDelta_dTau3 is not implemented for this backend"); };
232  virtual CoolPropDbl calc_d4alphar_dTau4(void){ throw NotImplementedError("calc_d4alphar_dTau4 is not implemented for this backend"); };
233 
234  // Derivatives of ideal-gas helmholtz energy
236  virtual CoolPropDbl calc_alpha0(void){ throw NotImplementedError("calc_alpha0 is not implemented for this backend"); };
238  virtual CoolPropDbl calc_dalpha0_dDelta(void){ throw NotImplementedError("calc_dalpha0_dDelta is not implemented for this backend"); };
240  virtual CoolPropDbl calc_dalpha0_dTau(void){ throw NotImplementedError("calc_dalpha0_dTau is not implemented for this backend"); };
242  virtual CoolPropDbl calc_d2alpha0_dDelta_dTau(void){ throw NotImplementedError("calc_d2alpha0_dDelta_dTau is not implemented for this backend"); };
244  virtual CoolPropDbl calc_d2alpha0_dDelta2(void){ throw NotImplementedError("calc_d2alpha0_dDelta2 is not implemented for this backend"); };
246  virtual CoolPropDbl calc_d2alpha0_dTau2(void){ throw NotImplementedError("calc_d2alpha0_dTau2 is not implemented for this backend"); };
248  virtual CoolPropDbl calc_d3alpha0_dDelta3(void){ throw NotImplementedError("calc_d3alpha0_dDelta3 is not implemented for this backend"); };
250  virtual CoolPropDbl calc_d3alpha0_dDelta2_dTau(void){ throw NotImplementedError("calc_d3alpha0_dDelta2_dTau is not implemented for this backend"); };
252  virtual CoolPropDbl calc_d3alpha0_dDelta_dTau2(void){ throw NotImplementedError("calc_d3alpha0_dDelta_dTau2 is not implemented for this backend"); };
254  virtual CoolPropDbl calc_d3alpha0_dTau3(void){ throw NotImplementedError("calc_d3alpha0_dTau3 is not implemented for this backend"); };
255 
256  virtual void calc_reducing_state(void){ throw NotImplementedError("calc_reducing_state is not implemented for this backend"); };
257 
259  virtual CoolPropDbl calc_Tmax(void){ throw NotImplementedError("calc_Tmax is not implemented for this backend"); };
261  virtual CoolPropDbl calc_Tmin(void){ throw NotImplementedError("calc_Tmin is not implemented for this backend"); };
263  virtual CoolPropDbl calc_pmax(void){ throw NotImplementedError("calc_pmax is not implemented for this backend"); };
264 
266  virtual CoolPropDbl calc_GWP20(void){ throw NotImplementedError("calc_GWP20 is not implemented for this backend"); };
268  virtual CoolPropDbl calc_GWP100(void){ throw NotImplementedError("calc_GWP100 is not implemented for this backend"); };
270  virtual CoolPropDbl calc_GWP500(void){ throw NotImplementedError("calc_GWP500 is not implemented for this backend"); };
272  virtual CoolPropDbl calc_ODP(void){ throw NotImplementedError("calc_ODP is not implemented for this backend"); };
274  virtual CoolPropDbl calc_flame_hazard(void){ throw NotImplementedError("calc_flame_hazard is not implemented for this backend"); };
276  virtual CoolPropDbl calc_health_hazard(void){ throw NotImplementedError("calc_health_hazard is not implemented for this backend"); };
278  virtual CoolPropDbl calc_physical_hazard(void){ throw NotImplementedError("calc_physical_hazard is not implemented for this backend"); };
280  virtual CoolPropDbl calc_dipole_moment(void){ throw NotImplementedError("calc_dipole_moment is not implemented for this backend"); };
281 
283  virtual CoolPropDbl calc_first_partial_deriv(parameters Of, parameters Wrt, parameters Constant);
285  virtual CoolPropDbl calc_second_partial_deriv(parameters Of1, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2);
286 
288  virtual CoolPropDbl calc_reduced_density(void){ throw NotImplementedError("calc_reduced_density is not implemented for this backend"); };
290  virtual CoolPropDbl calc_reciprocal_reduced_temperature(void){ throw NotImplementedError("calc_reciprocal_reduced_temperature is not implemented for this backend"); };
291 
293  virtual CoolPropDbl calc_Bvirial(void){ throw NotImplementedError("calc_Bvirial is not implemented for this backend"); };
295  virtual CoolPropDbl calc_Cvirial(void){ throw NotImplementedError("calc_Cvirial is not implemented for this backend"); };
297  virtual CoolPropDbl calc_dBvirial_dT(void){ throw NotImplementedError("calc_dBvirial_dT is not implemented for this backend"); };
299  virtual CoolPropDbl calc_dCvirial_dT(void){ throw NotImplementedError("calc_dCvirial_dT is not implemented for this backend"); };
301  virtual CoolPropDbl calc_compressibility_factor(void){ throw NotImplementedError("calc_compressibility_factor is not implemented for this backend"); };
302 
304  virtual std::string calc_name(void){ throw NotImplementedError("calc_name is not implemented for this backend"); };
306  virtual std::string calc_description(void){ throw NotImplementedError("calc_description is not implemented for this backend"); };
307 
309  virtual CoolPropDbl calc_Ttriple(void){ throw NotImplementedError("calc_Ttriple is not implemented for this backend"); };
311  virtual CoolPropDbl calc_p_triple(void){ throw NotImplementedError("calc_p_triple is not implemented for this backend"); };
312 
314  virtual CoolPropDbl calc_T_critical(void){ throw NotImplementedError("calc_T_critical is not implemented for this backend"); };
316  virtual CoolPropDbl calc_T_reducing(void){ throw NotImplementedError("calc_T_reducing is not implemented for this backend"); };
318  virtual CoolPropDbl calc_p_critical(void){ throw NotImplementedError("calc_p_critical is not implemented for this backend"); };
320  virtual CoolPropDbl calc_p_reducing(void){ throw NotImplementedError("calc_p_reducing is not implemented for this backend"); };
322  virtual CoolPropDbl calc_rhomolar_critical(void){ throw NotImplementedError("calc_rhomolar_critical is not implemented for this backend"); };
324  virtual CoolPropDbl calc_rhomass_critical(void){ throw NotImplementedError("calc_rhomass_critical is not implemented for this backend"); };
326  virtual CoolPropDbl calc_rhomolar_reducing(void){ throw NotImplementedError("calc_rhomolar_reducing is not implemented for this backend"); };
328  virtual void calc_phase_envelope(const std::string &type){ throw NotImplementedError("calc_phase_envelope is not implemented for this backend"); };
330  virtual CoolPropDbl calc_rhomass(void){ return rhomolar()*molar_mass(); }
331  virtual CoolPropDbl calc_hmass(void){ return hmolar() / molar_mass(); }
332  virtual CoolPropDbl calc_hmass_excess(void) { return hmolar_excess() / molar_mass(); }
333  virtual CoolPropDbl calc_smass(void){ return smolar() / molar_mass(); }
334  virtual CoolPropDbl calc_smass_excess(void) { return smolar_excess() / molar_mass(); }
335  virtual CoolPropDbl calc_cpmass(void){ return cpmolar() / molar_mass(); }
336  virtual CoolPropDbl calc_cp0mass(void){ return cp0molar() / molar_mass(); }
337  virtual CoolPropDbl calc_cvmass(void){ return cvmolar() / molar_mass(); }
338  virtual CoolPropDbl calc_umass(void){ return umolar() / molar_mass(); }
339  virtual CoolPropDbl calc_umass_excess(void) { return umolar_excess() / molar_mass(); }
340  virtual CoolPropDbl calc_gibbsmass(void){ return gibbsmolar() / molar_mass(); }
341  virtual CoolPropDbl calc_gibbsmass_excess(void) { return gibbsmolar_excess() / molar_mass(); }
342  virtual CoolPropDbl calc_helmholtzmass(void){ return helmholtzmolar() / molar_mass(); }
343  virtual CoolPropDbl calc_helmholtzmass_excess(void) { return helmholtzmolar_excess() / molar_mass(); }
344  virtual CoolPropDbl calc_volumemass_excess(void) { return volumemolar_excess() / molar_mass(); }
345 
347  virtual void update_states(void){ throw NotImplementedError("This backend does not implement update_states function"); };
348 
349  virtual CoolPropDbl calc_melting_line(int param, int given, CoolPropDbl value){ throw NotImplementedError("This backend does not implement calc_melting_line function"); };
350 
355  virtual CoolPropDbl calc_saturation_ancillary(parameters param, int Q, parameters given, double value){ throw NotImplementedError("This backend does not implement calc_saturation_ancillary"); };
356 
358  virtual phases calc_phase(void){ throw NotImplementedError("This backend does not implement calc_phase function"); };
360  virtual void calc_specify_phase(phases phase){ throw NotImplementedError("This backend does not implement calc_specify_phase function"); };
362  virtual void calc_unspecify_phase(void){ throw NotImplementedError("This backend does not implement calc_unspecify_phase function"); };
364  virtual std::vector<std::string> calc_fluid_names(void){ throw NotImplementedError("This backend does not implement calc_fluid_names function"); };
367  virtual const CoolProp::SimpleState & calc_state(const std::string &state){ throw NotImplementedError("calc_state is not implemented for this backend"); };
368 
369  virtual const CoolProp::PhaseEnvelopeData & calc_phase_envelope_data(void){ throw NotImplementedError("calc_phase_envelope_data is not implemented for this backend"); };
370 
371  virtual std::vector<CoolPropDbl> calc_mole_fractions_liquid(void){ throw NotImplementedError("calc_mole_fractions_liquid is not implemented for this backend"); };
372  virtual std::vector<CoolPropDbl> calc_mole_fractions_vapor(void){ throw NotImplementedError("calc_mole_fractions_vapor is not implemented for this backend"); };
373  virtual const std::vector<CoolPropDbl> calc_mass_fractions(void){ throw NotImplementedError("calc_mass_fractions is not implemented for this backend"); };
374 
376  virtual CoolPropDbl calc_fraction_min(void){ throw NotImplementedError("calc_fraction_min is not implemented for this backend"); };
378  virtual CoolPropDbl calc_fraction_max(void){ throw NotImplementedError("calc_fraction_max is not implemented for this backend"); };
379  virtual CoolPropDbl calc_T_freeze(void){ throw NotImplementedError("calc_T_freeze is not implemented for this backend"); };
380 
381  virtual CoolPropDbl calc_first_saturation_deriv(parameters Of1, parameters Wrt1){ throw NotImplementedError("calc_first_saturation_deriv is not implemented for this backend"); };
382  virtual CoolPropDbl calc_second_saturation_deriv(parameters Of1, parameters Wrt1, parameters Wrt2){ throw NotImplementedError("calc_second_saturation_deriv is not implemented for this backend"); };
383  virtual CoolPropDbl calc_first_two_phase_deriv(parameters Of, parameters Wrt, parameters Constant){ throw NotImplementedError("calc_first_two_phase_deriv is not implemented for this backend"); };
384  virtual CoolPropDbl calc_second_two_phase_deriv(parameters Of, parameters Wrt, parameters Constant, parameters Wrt2, parameters Constant2){ throw NotImplementedError("calc_second_two_phase_deriv is not implemented for this backend"); };
385  virtual CoolPropDbl calc_first_two_phase_deriv_splined(parameters Of, parameters Wrt, parameters Constant, CoolPropDbl x_end){ throw NotImplementedError("calc_first_two_phase_deriv_splined is not implemented for this backend"); };
386 
387  virtual CoolPropDbl calc_saturated_liquid_keyed_output(parameters key){ throw NotImplementedError("calc_saturated_liquid_keyed_output is not implemented for this backend"); };
388  virtual CoolPropDbl calc_saturated_vapor_keyed_output(parameters key){ throw NotImplementedError("calc_saturated_vapor_keyed_output is not implemented for this backend"); };
389  virtual void calc_ideal_curve(const std::string &type, std::vector<double> &T, std::vector<double> &p){ throw NotImplementedError("calc_ideal_curve is not implemented for this backend"); };
390 
392  virtual CoolPropDbl calc_T(void){ return _T; }
394  virtual CoolPropDbl calc_rhomolar(void){ return _rhomolar; }
395 
397  virtual double calc_tangent_plane_distance(const double T, const double p, const std::vector<double> &w, const double rhomolar_guess){ throw NotImplementedError("calc_tangent_plane_distance is not implemented for this backend"); };
398 
400  virtual void calc_true_critical_point(double &T, double &rho){ throw NotImplementedError("calc_true_critical_point is not implemented for this backend"); };
401 
402  virtual void calc_conformal_state(const std::string &reference_fluid, CoolPropDbl &T, CoolPropDbl &rhomolar){ throw NotImplementedError("calc_conformal_state is not implemented for this backend"); };
403 
404  virtual void calc_viscosity_contributions(CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical){ throw NotImplementedError("calc_viscosity_contributions is not implemented for this backend"); };
405  virtual void calc_conductivity_contributions(CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical){ throw NotImplementedError("calc_conductivity_contributions is not implemented for this backend"); };
406  virtual std::vector<CriticalState> calc_all_critical_points(void){ throw NotImplementedError("calc_all_critical_points is not implemented for this backend"); };
407  virtual void calc_build_spinodal(){ throw NotImplementedError("calc_build_spinodal is not implemented for this backend"); };
408  virtual SpinodalData calc_get_spinodal_data(){ throw NotImplementedError("calc_get_spinodal_data is not implemented for this backend"); };
409  virtual void calc_criticality_contour_values(double &L1star, double &M1star){ throw NotImplementedError("calc_criticality_contour_values is not implemented for this backend"); };
410 
412  virtual void mass_to_molar_inputs(CoolProp::input_pairs &input_pair, CoolPropDbl &value1, CoolPropDbl &value2);
413 
415  virtual void calc_change_EOS(const std::size_t i, const std::string &EOS_name){ throw NotImplementedError("calc_change_EOS is not implemented for this backend"); };
416 public:
417 
418  AbstractState() :_fluid_type(FLUID_TYPE_UNDEFINED), _phase(iphase_unknown){ clear(); }
419  virtual ~AbstractState(){};
420 
422 
428  static AbstractState * factory(const std::string &backend, const std::string &fluid_names)
429  {
430  return factory(backend, strsplit(fluid_names, '&'));
431  };
432 
456  static AbstractState * factory(const std::string &backend, const std::vector<std::string> &fluid_names);
457 
459  void set_T(CoolPropDbl T){ _T = T; }
460 
465  virtual std::string backend_name(void) = 0;
466 
467  // The derived classes must implement this function to define whether they use mole fractions (true) or mass fractions (false)
468  virtual bool using_mole_fractions(void) = 0;
469  virtual bool using_mass_fractions(void) = 0;
470  virtual bool using_volu_fractions(void) = 0;
471 
472  virtual void set_mole_fractions(const std::vector<CoolPropDbl> &mole_fractions) = 0;
473  virtual void set_mass_fractions(const std::vector<CoolPropDbl> &mass_fractions) = 0;
474  virtual void set_volu_fractions(const std::vector<CoolPropDbl> &mass_fractions){ throw NotImplementedError("Volume composition has not been implemented."); }
475 
476 
477 
497  virtual void set_reference_stateS(const std::string &reference_state){
498  throw NotImplementedError("Setting reference state has not been implemented for this backend. Try using CoolProp::set_reference_stateD instead.");
499  }
500 
506  virtual void set_reference_stateD(double T, double rhomolar, double hmolar0, double smolar0){
507  throw NotImplementedError("Setting reference state has not been implemented for this backend. Try using CoolProp::set_reference_stateD instead.");
508  }
509 
510 
511 #ifndef COOLPROPDBL_MAPS_TO_DOUBLE
512  void set_mole_fractions(const std::vector<double> &mole_fractions){ set_mole_fractions(std::vector<CoolPropDbl>(mole_fractions.begin(), mole_fractions.end())); };
513  void set_mass_fractions(const std::vector<double> &mass_fractions){ set_mass_fractions(std::vector<CoolPropDbl>(mass_fractions.begin(), mass_fractions.end())); };
514  void set_volu_fractions(const std::vector<double> &volu_fractions){ set_volu_fractions(std::vector<CoolPropDbl>(volu_fractions.begin(), volu_fractions.end())); };
515 #endif
516 
517  #ifdef EMSCRIPTEN
518  void set_mole_fractions_double(const std::vector<double> &mole_fractions){ set_mole_fractions(std::vector<CoolPropDbl>(mole_fractions.begin(), mole_fractions.end())); };
519  #endif
520 
522  std::vector<CoolPropDbl> mole_fractions_liquid(void){ return calc_mole_fractions_liquid(); };
524  std::vector<double> mole_fractions_liquid_double(void){
525  std::vector<CoolPropDbl> x = calc_mole_fractions_liquid();
526  return std::vector<double>(x.begin(), x.end());
527  };
528 
530  std::vector<CoolPropDbl> mole_fractions_vapor(void){ return calc_mole_fractions_vapor(); };
532  std::vector<double> mole_fractions_vapor_double(void){
533  std::vector<CoolPropDbl> y = calc_mole_fractions_vapor();
534  return std::vector<double>(y.begin(), y.end());
535  };
536 
538  virtual const std::vector<CoolPropDbl> & get_mole_fractions(void) = 0;
540  virtual const std::vector<CoolPropDbl> get_mass_fractions(void){
541  return this->calc_mass_fractions();
542  };
543 
545  virtual void update(CoolProp::input_pairs input_pair, double Value1, double Value2) = 0;
546 
549  virtual void update_with_guesses(CoolProp::input_pairs input_pair, double Value1, double Value2, const GuessesStructure &guesses){ throw NotImplementedError("update_with_guesses is not implemented for this backend"); };
550 
554  virtual bool available_in_high_level(void){ return true; }
555 
557  virtual std::string fluid_param_string(const std::string &){ throw NotImplementedError("fluid_param_string has not been implemented for this backend"); }
558 
560  std::vector<std::string> fluid_names(void);
561 
566  virtual const double get_fluid_constant(std::size_t i, parameters param) const{ throw NotImplementedError("get_fluid_constant is not implemented for this backend"); };
567 ;
568 
570  virtual void set_binary_interaction_double(const std::string &CAS1, const std::string &CAS2, const std::string &parameter, const double value){ throw NotImplementedError("set_binary_interaction_double is not implemented for this backend"); };
572  virtual void set_binary_interaction_double(const std::size_t i, const std::size_t j, const std::string &parameter, const double value){ throw NotImplementedError("set_binary_interaction_double is not implemented for this backend"); };
574  virtual void set_binary_interaction_string(const std::string &CAS1, const std::string &CAS2, const std::string &parameter, const std::string &value){ throw NotImplementedError("set_binary_interaction_string is not implemented for this backend"); };
576  virtual void set_binary_interaction_string(const std::size_t i, const std::size_t j, const std::string &parameter, const std::string &value){ throw NotImplementedError("set_binary_interaction_string is not implemented for this backend"); };
578  virtual double get_binary_interaction_double(const std::string &CAS1, const std::string &CAS2, const std::string &parameter){ throw NotImplementedError("get_binary_interaction_double is not implemented for this backend"); };
580  virtual double get_binary_interaction_double(const std::size_t i, const std::size_t j, const std::string &parameter){ throw NotImplementedError("get_binary_interaction_double is not implemented for this backend"); };
582  virtual std::string get_binary_interaction_string(const std::string &CAS1, const std::string &CAS2, const std::string &parameter){ throw NotImplementedError("get_binary_interaction_string is not implemented for this backend"); };
584  virtual void apply_simple_mixing_rule(std::size_t i, std::size_t j, const std::string &model) { throw NotImplementedError("apply_simple_mixing_rule is not implemented for this backend"); };
586  virtual void set_cubic_alpha_C(const size_t i, const std::string &parameter, const double c1, const double c2, const double c3) { throw ValueError("set_cubic_alpha_C only defined for cubic backends"); };
588  virtual void set_fluid_parameter_double(const size_t i, const std::string &parameter, const double value) { throw ValueError("set_fluid_parameter_double only defined for cubic backends"); };
590  virtual double get_fluid_parameter_double(const size_t i, const std::string &parameter) { throw ValueError("get_fluid_parameter_double only defined for cubic backends"); };
591 
593  virtual bool clear();
595  virtual bool clear_comp_change();
596 
597 
601  virtual const CoolProp::SimpleState & get_reducing_state(){ return _reducing; };
602 
604  const CoolProp::SimpleState & get_state(const std::string &state){ return calc_state(state); };
605 
607  double Tmin(void);
609  double Tmax(void);
611  double pmax(void);
613  double Ttriple(void);
614 
616  phases phase(void){ return calc_phase(); };
618  void specify_phase(phases phase){ calc_specify_phase(phase); };
620  void unspecify_phase(void){ calc_unspecify_phase(); };
621 
623  double T_critical(void);
625  double p_critical(void);
627  double rhomolar_critical(void);
629  double rhomass_critical(void);
630 
632  std::vector<CriticalState> all_critical_points(void){ return calc_all_critical_points(); };
633 
635  void build_spinodal(){ calc_build_spinodal(); };
636 
638  SpinodalData get_spinodal_data(){ return calc_get_spinodal_data(); };
639 
641  void criticality_contour_values(double &L1star, double &M1star){ return calc_criticality_contour_values(L1star, M1star); }
642 
666  double tangent_plane_distance(const double T, const double p, const std::vector<double> &w, const double rhomolar_guess = -1){ return calc_tangent_plane_distance(T, p, w, rhomolar_guess); };
667 
669  double T_reducing(void);
671  double rhomolar_reducing(void);
673  double rhomass_reducing(void);
674 
676  double p_triple(void);
677 
679  std::string name(){ return calc_name(); };
681  std::string description(){ return calc_description(); };
682 
684  double dipole_moment(){ return calc_dipole_moment(); }
685 
686  // ----------------------------------------
687  // Bulk properties - temperature and density are directly calculated every time
688  // All other parameters are calculated on an as-needed basis
689  // ----------------------------------------
691  double keyed_output(parameters key);
693  double trivial_keyed_output(parameters key);
695  double saturated_liquid_keyed_output(parameters key){ return calc_saturated_liquid_keyed_output(key); };
697  double saturated_vapor_keyed_output(parameters key){ return calc_saturated_vapor_keyed_output(key); };
698 
700  double T(void) { return calc_T(); };
702  double rhomolar(void){ return calc_rhomolar(); };
704  double rhomass(void){ return calc_rhomass(); };
706  double p(void) { return _p; };
708  double Q(void) { return _Q; };
710  double tau(void);
712  double delta(void);
714  double molar_mass(void);
716  double acentric_factor(void);
718  double gas_constant(void);
720  double Bvirial(void);
722  double dBvirial_dT(void);
724  double Cvirial(void);
726  double dCvirial_dT(void);
728  double compressibility_factor(void);
730  double hmolar(void);
732  double hmolar_residual(void);
734  double hmass(void){ return calc_hmass(); };
736  double hmolar_excess(void);
738  double hmass_excess(void) { return calc_hmass_excess(); };
740  double smolar(void);
742  double smolar_residual(void);
744  double smass(void){ return calc_smass(); };
746  double smolar_excess(void);
748  double smass_excess(void) { return calc_smass_excess(); };
750  double umolar(void);
752  double umass(void){ return calc_umass(); };
754  double umolar_excess(void);
756  double umass_excess(void) { return calc_umass_excess(); };
758  double cpmolar(void);
760  double cpmass(void){ return calc_cpmass(); };
762  double cp0molar(void);
764  double cp0mass(void){ return calc_cp0mass(); };
766  double cvmolar(void);
768  double cvmass(void){ return calc_cvmass(); };
770  double gibbsmolar(void);
772  double gibbsmolar_residual(void);
774  double gibbsmass(void){ return calc_gibbsmass(); };
776  double gibbsmolar_excess(void);
778  double gibbsmass_excess(void) { return calc_gibbsmass_excess(); };
780  double helmholtzmolar(void);
782  double helmholtzmass(void){ return calc_helmholtzmass(); };
784  double helmholtzmolar_excess(void);
786  double helmholtzmass_excess(void) { return calc_helmholtzmass_excess(); };
788  double volumemolar_excess(void);
790  double volumemass_excess(void) { return calc_volumemass_excess(); };
792  double speed_sound(void);
794  double isothermal_compressibility(void);
796  double isobaric_expansion_coefficient(void);
798  double isentropic_expansion_coefficient(void);
800  double fugacity_coefficient(std::size_t i);
802  std::vector<double> fugacity_coefficients();
804  double fugacity(std::size_t i);
806  double chemical_potential(std::size_t i);
815  double fundamental_derivative_of_gas_dynamics(void);
817  double PIP(){ return calc_PIP(); };
818 
820  void true_critical_point(double &T, double &rho){ calc_true_critical_point(T, rho); }
821 
829  void ideal_curve(const std::string &type, std::vector<double> &T, std::vector<double> &p){ calc_ideal_curve(type, T, p); };
830 
831  // ----------------------------------------
832  // Partial derivatives
833  // ----------------------------------------
834 
839  CoolPropDbl first_partial_deriv(parameters Of, parameters Wrt, parameters Constant){return calc_first_partial_deriv(Of, Wrt, Constant);};
840 
864  CoolPropDbl second_partial_deriv(parameters Of1, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2){return calc_second_partial_deriv(Of1,Wrt1,Constant1,Wrt2,Constant2);};
865 
887  CoolPropDbl first_saturation_deriv(parameters Of1, parameters Wrt1){return calc_first_saturation_deriv(Of1,Wrt1);};
888 
906  CoolPropDbl second_saturation_deriv(parameters Of1, parameters Wrt1, parameters Wrt2){return calc_second_saturation_deriv(Of1,Wrt1,Wrt2);};
907 
928  return calc_first_two_phase_deriv(Of, Wrt, Constant);
929  };
930 
947  double second_two_phase_deriv(parameters Of, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2){
948  return calc_second_two_phase_deriv(Of, Wrt1, Constant1, Wrt2, Constant2);
949  };
950 
971  double first_two_phase_deriv_splined(parameters Of, parameters Wrt, parameters Constant, double x_end){
972  return calc_first_two_phase_deriv_splined(Of, Wrt, Constant, x_end);
973  };
974 
975  // ----------------------------------------
976  // Phase envelope for mixtures
977  // ----------------------------------------
978 
984  void build_phase_envelope(const std::string &type = "");
988  const CoolProp::PhaseEnvelopeData &get_phase_envelope_data(){return calc_phase_envelope_data();};
989 
990  // ----------------------------------------
991  // Ancillary equations
992  // ----------------------------------------
993 
995  virtual bool has_melting_line(void){return false;};
1000  double melting_line(int param, int given, double value);
1006  double saturation_ancillary(parameters param, int Q, parameters given, double value);
1007 
1008  // ----------------------------------------
1009  // Transport properties
1010  // ----------------------------------------
1012  double viscosity(void);
1014  void viscosity_contributions(CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical){ calc_viscosity_contributions(dilute, initial_density, residual, critical); };
1016  double conductivity(void);
1018  void conductivity_contributions(CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical){ calc_conductivity_contributions(dilute, initial_density, residual, critical); };
1020  double surface_tension(void);
1022  double Prandtl(void){return cpmass()*viscosity()/conductivity();};
1029  void conformal_state(const std::string &reference_fluid, CoolPropDbl &T, CoolPropDbl &rhomolar){
1030  return calc_conformal_state(reference_fluid, T, rhomolar);
1031  };
1032 
1037  void change_EOS(const std::size_t i, const std::string &EOS_name){ calc_change_EOS(i, EOS_name); }
1038 
1039  // ----------------------------------------
1040  // Helmholtz energy and derivatives
1041  // ----------------------------------------
1043  CoolPropDbl alpha0(void){
1044  if (!_alpha0) _alpha0 = calc_alpha0();
1045  return _alpha0;
1046  };
1048  CoolPropDbl dalpha0_dDelta(void){
1049  if (!_dalpha0_dDelta) _dalpha0_dDelta = calc_dalpha0_dDelta();
1050  return _dalpha0_dDelta;
1051  };
1053  CoolPropDbl dalpha0_dTau(void){
1054  if (!_dalpha0_dTau) _dalpha0_dTau = calc_dalpha0_dTau();
1055  return _dalpha0_dTau;
1056  };
1058  CoolPropDbl d2alpha0_dDelta2(void){
1059  if (!_d2alpha0_dDelta2) _d2alpha0_dDelta2 = calc_d2alpha0_dDelta2();
1060  return _d2alpha0_dDelta2;
1061  };
1063  CoolPropDbl d2alpha0_dDelta_dTau(void){
1064  if (!_d2alpha0_dDelta_dTau) _d2alpha0_dDelta_dTau = calc_d2alpha0_dDelta_dTau();
1065  return _d2alpha0_dDelta_dTau;
1066  };
1068  CoolPropDbl d2alpha0_dTau2(void){
1069  if (!_d2alpha0_dTau2) _d2alpha0_dTau2 = calc_d2alpha0_dTau2();
1070  return _d2alpha0_dTau2;
1071  };
1073  CoolPropDbl d3alpha0_dTau3(void){
1074  if (!_d3alpha0_dTau3) _d3alpha0_dTau3 = calc_d3alpha0_dTau3();
1075  return _d3alpha0_dTau3;
1076  };
1078  CoolPropDbl d3alpha0_dDelta_dTau2(void){
1079  if (!_d3alpha0_dDelta_dTau2) _d3alpha0_dDelta_dTau2 = calc_d3alpha0_dDelta_dTau2();
1080  return _d3alpha0_dDelta_dTau2;
1081  };
1083  CoolPropDbl d3alpha0_dDelta2_dTau(void){
1084  if (!_d3alpha0_dDelta2_dTau) _d3alpha0_dDelta2_dTau = calc_d3alpha0_dDelta2_dTau();
1085  return _d3alpha0_dDelta2_dTau;
1086  };
1088  CoolPropDbl d3alpha0_dDelta3(void){
1089  if (!_d3alpha0_dDelta3) _d3alpha0_dDelta3 = calc_d3alpha0_dDelta3();
1090  return _d3alpha0_dDelta3;
1091  };
1092 
1094  CoolPropDbl alphar(void){
1095  if (!_alphar) _alphar = calc_alphar();
1096  return _alphar;
1097  };
1099  CoolPropDbl dalphar_dDelta(void){
1100  if (!_dalphar_dDelta) _dalphar_dDelta = calc_dalphar_dDelta();
1101  return _dalphar_dDelta;
1102  };
1104  CoolPropDbl dalphar_dTau(void){
1105  if (!_dalphar_dTau) _dalphar_dTau = calc_dalphar_dTau();
1106  return _dalphar_dTau;
1107  };
1109  CoolPropDbl d2alphar_dDelta2(void){
1110  if (!_d2alphar_dDelta2) _d2alphar_dDelta2 = calc_d2alphar_dDelta2();
1111  return _d2alphar_dDelta2;
1112  };
1114  CoolPropDbl d2alphar_dDelta_dTau(void){
1115  if (!_d2alphar_dDelta_dTau) _d2alphar_dDelta_dTau = calc_d2alphar_dDelta_dTau();
1116  return _d2alphar_dDelta_dTau;
1117  };
1119  CoolPropDbl d2alphar_dTau2(void){
1120  if (!_d2alphar_dTau2) _d2alphar_dTau2 = calc_d2alphar_dTau2();
1121  return _d2alphar_dTau2;
1122  };
1124  CoolPropDbl d3alphar_dDelta3(void){
1125  if (!_d3alphar_dDelta3) _d3alphar_dDelta3 = calc_d3alphar_dDelta3();
1126  return _d3alphar_dDelta3;
1127  };
1129  CoolPropDbl d3alphar_dDelta2_dTau(void){
1130  if (!_d3alphar_dDelta2_dTau) _d3alphar_dDelta2_dTau = calc_d3alphar_dDelta2_dTau();
1131  return _d3alphar_dDelta2_dTau;
1132  };
1134  CoolPropDbl d3alphar_dDelta_dTau2(void){
1135  if (!_d3alphar_dDelta_dTau2) _d3alphar_dDelta_dTau2 = calc_d3alphar_dDelta_dTau2();
1136  return _d3alphar_dDelta_dTau2;
1137  };
1139  CoolPropDbl d3alphar_dTau3(void){
1140  if (!_d3alphar_dTau3) _d3alphar_dTau3 = calc_d3alphar_dTau3();
1141  return _d3alphar_dTau3;
1142  };
1144  CoolPropDbl d4alphar_dDelta4(void){
1145  if (!_d4alphar_dDelta4) _d4alphar_dDelta4 = calc_d4alphar_dDelta4();
1146  return _d4alphar_dDelta4;
1147  };
1149  CoolPropDbl d4alphar_dDelta3_dTau(void){
1150  if (!_d4alphar_dDelta3_dTau) _d4alphar_dDelta3_dTau = calc_d4alphar_dDelta3_dTau();
1151  return _d4alphar_dDelta3_dTau;
1152  };
1154  CoolPropDbl d4alphar_dDelta2_dTau2(void){
1155  if (!_d4alphar_dDelta2_dTau2) _d4alphar_dDelta2_dTau2 = calc_d4alphar_dDelta2_dTau2();
1156  return _d4alphar_dDelta2_dTau2;
1157  };
1159  CoolPropDbl d4alphar_dDelta_dTau3(void){
1160  if (!_d4alphar_dDelta_dTau3) _d4alphar_dDelta_dTau3 = calc_d4alphar_dDelta_dTau3();
1161  return _d4alphar_dDelta_dTau3;
1162  };
1164  CoolPropDbl d4alphar_dTau4(void){
1165  if (!_d4alphar_dTau4) _d4alphar_dTau4 = calc_d4alphar_dTau4();
1166  return _d4alphar_dTau4;
1167  };
1168 };
1169 
1178 public:
1179  virtual AbstractState * get_AbstractState(const std::vector<std::string> &fluid_names) = 0;
1180  virtual ~AbstractStateGenerator() {};
1181 };
1182 
1186 void register_backend(const backend_families &bf, shared_ptr<AbstractStateGenerator> gen);
1187 
1188 template <class T>
1190 public:
1192  register_backend(bf, shared_ptr<AbstractStateGenerator>(new T()));
1193  };
1194 };
1195 
1196 
1197 } /* namespace CoolProp */
1198 #endif /* ABSTRACTSTATE_H_ */
double T(void)
Return the temperature in K.
Definition: AbstractState.h:700
phases imposed_phase_index
If the phase is imposed, the imposed phase index.
Definition: AbstractState.h:76
virtual CoolPropDbl calc_dalphar_dDelta(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:205
std::vector< double > delta
The reduced density ( )
Definition: AbstractState.h:25
CoolPropDbl dalphar_dDelta(void)
Return the term .
Definition: AbstractState.h:1099
virtual std::string get_binary_interaction_string(const std::string &CAS1, const std::string &CAS2, const std::string &parameter)
Get binary mixture string value (EXPERT USE ONLY!!!)
Definition: AbstractState.h:582
double cvmass(void)
Return the mass constant volume specific heat in J/kg/K.
Definition: AbstractState.h:768
virtual CoolPropDbl calc_hmolar_residual(void)
Using this backend, calculate the residual molar enthalpy in J/mol.
Definition: AbstractState.h:145
CoolPropDbl d4alphar_dDelta2_dTau2(void)
Return the term .
Definition: AbstractState.h:1154
Supercritical liquid (p > pc, T < Tc)
Definition: DataStructures.h:161
virtual CoolPropDbl calc_physical_hazard(void)
Using this backend, calculate the physical hazard.
Definition: AbstractState.h:278
SimpleState _critical
Two important points.
Definition: AbstractState.h:91
double T
temperature in K
Definition: AbstractState.h:34
static AbstractState * factory(const std::string &backend, const std::string &fluid_names)
A factory function to return a pointer to a new-allocated instance of one of the backends.
Definition: AbstractState.h:428
double rhomolar(void)
Return the molar density in mol/m^3.
Definition: AbstractState.h:702
virtual CoolPropDbl calc_flame_hazard(void)
Using this backend, calculate the flame hazard.
Definition: AbstractState.h:274
CoolPropDbl first_partial_deriv(parameters Of, parameters Wrt, parameters Constant)
The first partial derivative in homogeneous phases.
Definition: AbstractState.h:839
virtual void apply_simple_mixing_rule(std::size_t i, std::size_t j, const std::string &model)
Apply a simple mixing rule (EXPERT USE ONLY!!!)
Definition: AbstractState.h:584
virtual CoolPropDbl calc_ODP(void)
Using this backend, calculate the ozone depletion potential (ODP)
Definition: AbstractState.h:272
std::vector< double > tau
The reciprocal reduced temperature ( )
Definition: AbstractState.h:25
CoolPropDbl dalpha0_dTau(void)
Return the term .
Definition: AbstractState.h:1053
CoolPropDbl d2alpha0_dDelta2(void)
Return the term .
Definition: AbstractState.h:1058
std::string name()
Return the name - backend dependent.
Definition: AbstractState.h:679
virtual CoolPropDbl calc_Bvirial(void)
Using this backend, calculate the second virial coefficient.
Definition: AbstractState.h:293
virtual CoolPropDbl calc_dalpha0_dDelta(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:238
double tangent_plane_distance(const double T, const double p, const std::vector< double > &w, const double rhomolar_guess=-1)
Return the tangent plane distance for a given trial composition w.
Definition: AbstractState.h:666
virtual CoolPropDbl calc_saturation_ancillary(parameters param, int Q, parameters given, double value)
Definition: AbstractState.h:355
virtual CoolPropDbl calc_d2alpha0_dTau2(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:246
CoolPropDbl d3alpha0_dDelta3(void)
Return the term .
Definition: AbstractState.h:1088
virtual CoolPropDbl calc_d2alphar_dTau2(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:213
virtual CoolPropDbl calc_reduced_density(void)
Using this backend, calculate the reduced density (rho/rhoc)
Definition: AbstractState.h:288
virtual CoolPropDbl calc_hmolar(void)
Using this backend, calculate the molar enthalpy in J/mol.
Definition: AbstractState.h:143
double hmass(void)
Return the mass enthalpy in J/kg.
Definition: AbstractState.h:734
virtual CoolPropDbl calc_d3alphar_dDelta_dTau2(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:219
CachedElement _viscosity
Transport properties.
Definition: AbstractState.h:105
virtual void update_with_guesses(CoolProp::input_pairs input_pair, double Value1, double Value2, const GuessesStructure &guesses)
Update the state using two state variables and providing guess values Some or all of the guesses will...
Definition: AbstractState.h:549
virtual CoolPropDbl calc_cpmolar_idealgas(void)
Using this backend, calculate the ideal gas molar constant-pressure specific heat in J/mol/K...
Definition: AbstractState.h:155
double gibbsmass(void)
Return the Gibbs energy in J/kg.
Definition: AbstractState.h:774
Subcritical liquid.
Definition: DataStructures.h:158
virtual CoolPropDbl calc_d3alphar_dTau3(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:221
virtual void calc_specify_phase(phases phase)
Using this backend, specify the phase to be used for all further calculations.
Definition: AbstractState.h:360
virtual CoolPropDbl calc_acentric_factor(void)
Using this backend, calculate the acentric factor.
Definition: AbstractState.h:181
Supercritical gas (p < pc, T > Tc)
Definition: DataStructures.h:160
virtual std::vector< std::string > calc_fluid_names(void)
Using this backend, get a vector of fluid names.
Definition: AbstractState.h:364
phases
These are constants for the phases of the fluid.
Definition: DataStructures.h:158
double rhomass(void)
Return the mass density in kg/m^3.
Definition: AbstractState.h:704
double smass_excess(void)
Return the molar entropy in J/kg/K.
Definition: AbstractState.h:748
phases phase(void)
Get the phase of the state.
Definition: AbstractState.h:616
virtual CoolPropDbl calc_compressibility_factor(void)
Using this backend, calculate the compressibility factor Z .
Definition: AbstractState.h:301
CoolPropDbl d3alpha0_dTau3(void)
Return the term .
Definition: AbstractState.h:1073
virtual void calc_excess_properties(void)
Using this backend, calculate and cache the excess properties.
Definition: AbstractState.h:199
virtual std::string calc_description(void)
Using this backend, get the description of the fluid.
Definition: AbstractState.h:306
virtual CoolPropDbl calc_dalpha0_dTau(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:240
double cp0mass(void)
Return the mass constant pressure specific heat for ideal gas part only in J/kg/K.
Definition: AbstractState.h:764
virtual void set_binary_interaction_string(const std::size_t i, const std::size_t j, const std::string &parameter, const std::string &value)
Set binary mixture string parameter (EXPERT USE ONLY!!!)
Definition: AbstractState.h:576
virtual CoolPropDbl calc_fraction_min(void)
Get the minimum fraction (mole, mass, volume) for incompressible fluid.
Definition: AbstractState.h:376
virtual CoolPropDbl calc_d4alphar_dDelta2_dTau2(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:228
double hmass_excess(void)
Return the excess mass enthalpy in J/kg.
Definition: AbstractState.h:738
long _fluid_type
Some administrative variables.
Definition: AbstractState.h:74
virtual bool has_melting_line(void)
Return true if the fluid has a melting line - default is false, but can be re-implemented by derived ...
Definition: AbstractState.h:995
void conformal_state(const std::string &reference_fluid, CoolPropDbl &T, CoolPropDbl &rhomolar)
Find the conformal state needed for ECS.
Definition: AbstractState.h:1029
virtual CoolPropDbl calc_health_hazard(void)
Using this backend, calculate the health hazard.
Definition: AbstractState.h:276
virtual CoolPropDbl calc_Cvirial(void)
Using this backend, calculate the third virial coefficient.
Definition: AbstractState.h:295
Unknown phase.
Definition: DataStructures.h:165
std::vector< CriticalState > all_critical_points(void)
Return the vector of critical points, including points that are unstable or correspond to negative pr...
Definition: AbstractState.h:632
An abstract AbstractState generator class.
Definition: AbstractState.h:1177
double cpmass(void)
Return the mass constant pressure specific heat in J/kg/K.
Definition: AbstractState.h:760
CoolPropDbl d2alphar_dDelta_dTau(void)
Return the term .
Definition: AbstractState.h:1114
virtual CoolPropDbl calc_Tmin(void)
Using this backend, calculate the minimum temperature in K.
Definition: AbstractState.h:261
double umass_excess(void)
Return the excess internal energy in J/kg.
Definition: AbstractState.h:756
virtual CoolPropDbl calc_T_reducing(void)
Using this backend, get the reducing point temperature in K.
Definition: AbstractState.h:316
The mother of all state classes.
Definition: AbstractState.h:70
virtual CoolPropDbl calc_reciprocal_reduced_temperature(void)
Using this backend, calculate the reciprocal reduced temperature (Tc/T)
Definition: AbstractState.h:290
virtual CoolPropDbl calc_d2alpha0_dDelta_dTau(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:242
CoolPropDbl d4alphar_dDelta3_dTau(void)
Return the term .
Definition: AbstractState.h:1149
virtual CoolPropDbl calc_GWP500(void)
Using this backend, calculate the 500-year global warming potential (GWP)
Definition: AbstractState.h:270
virtual std::string fluid_param_string(const std::string &)
Return a string from the backend for the mixture/fluid - backend dependent - could be CAS #...
Definition: AbstractState.h:557
CoolPropDbl d3alpha0_dDelta_dTau2(void)
Return the term .
Definition: AbstractState.h:1078
virtual void set_binary_interaction_double(const std::string &CAS1, const std::string &CAS2, const std::string &parameter, const double value)
Set binary mixture floating point parameter (EXPERT USE ONLY!!!)
Definition: AbstractState.h:570
virtual const double get_fluid_constant(std::size_t i, parameters param) const
Get a constant for one of the fluids forming this mixture.
Definition: AbstractState.h:566
virtual CoolPropDbl calc_d2alphar_dDelta_dTau(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:211
double helmholtzmass(void)
Return the Helmholtz energy in J/kg.
Definition: AbstractState.h:782
std::vector< CoolPropDbl > mole_fractions_liquid(void)
Get the mole fractions of the equilibrium liquid phase.
Definition: AbstractState.h:522
void set_T(CoolPropDbl T)
Set the internal variable T without a flash call (expert use only!)
Definition: AbstractState.h:459
virtual void calc_true_critical_point(double &T, double &rho)
Using this backend, return true critical point where dp/drho|T = 0 and d2p/drho^2|T = 0...
Definition: AbstractState.h:400
CoolPropDbl d3alphar_dTau3(void)
Return the term .
Definition: AbstractState.h:1139
CachedElement _rho_spline
Smoothing values.
Definition: AbstractState.h:121
double first_two_phase_deriv_splined(parameters Of, parameters Wrt, parameters Constant, double x_end)
Calculate the first "two-phase" derivative as described by Thorade and Sadaat, EAS, 2013.
Definition: AbstractState.h:971
virtual CoolPropDbl calc_pmax(void)
Using this backend, calculate the maximum pressure in Pa.
Definition: AbstractState.h:263
This structure holds values obtained while tracing the spinodal curve (most often in the process of f...
Definition: AbstractState.h:23
A data structure to hold the data for a phase envelope.
Definition: PhaseEnvelope.h:14
virtual CoolPropDbl calc_chemical_potential(std::size_t i)
Using this backend, calculate the chemical potential in J/mol.
Definition: AbstractState.h:193
virtual CoolPropDbl calc_Ttriple(void)
Using this backend, get the triple point temperature in K.
Definition: AbstractState.h:309
backend_families
The structure is taken directly from the AbstractState class.
Definition: DataStructures.h:409
At the critical point.
Definition: DataStructures.h:162
void true_critical_point(double &T, double &rho)
Calculate the "true" critical point for pure fluids where dpdrho|T and d2p/drho2|T are equal to zero...
Definition: AbstractState.h:820
virtual void set_binary_interaction_string(const std::string &CAS1, const std::string &CAS2, const std::string &parameter, const std::string &value)
Set binary mixture string parameter (EXPERT USE ONLY!!!)
Definition: AbstractState.h:574
double saturated_liquid_keyed_output(parameters key)
Get an output from the saturated liquid state by key.
Definition: AbstractState.h:695
CoolPropDbl first_saturation_deriv(parameters Of1, parameters Wrt1)
The first partial derivative along the saturation curve.
Definition: AbstractState.h:887
std::string description()
Return the description - backend dependent.
Definition: AbstractState.h:681
virtual CoolPropDbl calc_viscosity(void)
Using this backend, calculate the viscosity in Pa-s.
Definition: AbstractState.h:173
void conductivity_contributions(CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical)
Return the thermal conductivity contributions, each in W/m/K.
Definition: AbstractState.h:1018
void specify_phase(phases phase)
Specify the phase for all further calculations with this state class.
Definition: AbstractState.h:618
virtual void calc_change_EOS(const std::size_t i, const std::string &EOS_name)
Change the equation of state for a given component to a specified EOS.
Definition: AbstractState.h:415
CachedElement _alpha0
Cached low-level elements for in-place calculation of other properties.
Definition: AbstractState.h:125
input_pairs
These are input pairs that can be used for the update function (in each pair, input keys are sorted a...
Definition: DataStructures.h:216
virtual CoolPropDbl calc_d3alpha0_dTau3(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:254
virtual void set_cubic_alpha_C(const size_t i, const std::string &parameter, const double c1, const double c2, const double c3)
Set the cubic alpha function&#39;s constants:
Definition: AbstractState.h:586
void build_spinodal()
Construct the spinodal curve for the mixture (or pure fluid)
Definition: AbstractState.h:635
virtual void calc_phase_envelope(const std::string &type)
Using this backend, construct the phase envelope, the variable type describes the type of phase envel...
Definition: AbstractState.h:328
virtual bool available_in_high_level(void)
A function that says whether the backend instance can be instantiated in the high-level interface In ...
Definition: AbstractState.h:554
virtual void calc_unspecify_phase(void)
Using this backend, unspecify the phase.
Definition: AbstractState.h:362
Subcritical gas.
Definition: DataStructures.h:163
Definition: Exceptions.h:26
CachedElement _hmolar_excess
Excess properties.
Definition: AbstractState.h:113
CoolPropDbl d2alpha0_dTau2(void)
Return the term .
Definition: AbstractState.h:1068
virtual CoolPropDbl calc_Tmax(void)
Using this backend, calculate the maximum temperature in K.
Definition: AbstractState.h:259
virtual CoolPropDbl calc_conductivity(void)
Using this backend, calculate the thermal conductivity in W/m/K.
Definition: AbstractState.h:175
void change_EOS(const std::size_t i, const std::string &EOS_name)
Change the equation of state for a given component to a specified EOS.
Definition: AbstractState.h:1037
CoolPropDbl d3alphar_dDelta3(void)
Return the term .
Definition: AbstractState.h:1124
double helmholtzmass_excess(void)
Return the excess Helmholtz energy in J/kg.
Definition: AbstractState.h:786
CoolPropDbl d3alphar_dDelta2_dTau(void)
Return the term .
Definition: AbstractState.h:1129
virtual CoolPropDbl calc_isobaric_expansion_coefficient(void)
Using this backend, calculate the isobaric expansion coefficient in 1/K.
Definition: AbstractState.h:169
double PIP()
Return the phase identification parameter (PIP) of G. Venkatarathnam and L.R. Oellrich, "Identification of the phase of a fluid using partial derivatives of pressure, volume, and temperature without reference to saturation properties: Applications in phase equilibria calculations".
Definition: AbstractState.h:817
virtual CoolPropDbl calc_molar_mass(void)
Using this backend, calculate the molar mass in kg/mol.
Definition: AbstractState.h:179
double Q(void)
Return the vapor quality (mol/mol); Q = 0 for saturated liquid.
Definition: AbstractState.h:708
CachedElement _rhoLmolar
Two-Phase variables.
Definition: AbstractState.h:136
virtual CoolPropDbl calc_dBvirial_dT(void)
Using this backend, calculate the derivative dB/dT.
Definition: AbstractState.h:297
virtual CoolPropDbl calc_umolar(void)
Using this backend, calculate the molar internal energy in J/mol.
Definition: AbstractState.h:151
virtual CoolPropDbl calc_T(void)
Using this backend, get the temperature.
Definition: AbstractState.h:392
virtual double calc_tangent_plane_distance(const double T, const double p, const std::vector< double > &w, const double rhomolar_guess)
Using this backend, calculate the tangent plane distance for a given trial composition.
Definition: AbstractState.h:397
CoolPropDbl d4alphar_dDelta_dTau3(void)
Return the term .
Definition: AbstractState.h:1159
virtual CoolPropDbl calc_d3alpha0_dDelta3(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:248
virtual CoolPropDbl calc_helmholtzmolar(void)
Using this backend, calculate the molar Helmholtz energy in J/mol.
Definition: AbstractState.h:163
CachedElement _molar_mass
Molar mass [mol/kg].
Definition: AbstractState.h:94
double saturated_vapor_keyed_output(parameters key)
Get an output from the saturated vapor state by key.
Definition: AbstractState.h:697
std::vector< double > mole_fractions_vapor_double(void)
Get the mole fractions of the equilibrium vapor phase (but as a double for use in SWIG wrapper) ...
Definition: AbstractState.h:532
virtual CoolPropDbl calc_gas_constant(void)
Using this backend, calculate the universal gas constant in J/mol/K.
Definition: AbstractState.h:185
virtual void update_states(void)
Update the states after having changed the reference state for enthalpy and entropy.
Definition: AbstractState.h:347
virtual CoolPropDbl calc_d2alpha0_dDelta2(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:244
void criticality_contour_values(double &L1star, double &M1star)
Calculate the criticality contour values and .
Definition: AbstractState.h:641
virtual double get_fluid_parameter_double(const size_t i, const std::string &parameter)
Double fluid parameter (currently the volume translation parameter for cubic)
Definition: AbstractState.h:590
virtual CoolPropDbl calc_isentropic_expansion_coefficient(void)
Using this backend, calculate the isentropic expansion coefficient .
Definition: AbstractState.h:171
CachedElement _gas_constant
Universal gas constant [J/mol/K].
Definition: AbstractState.h:97
virtual CoolPropDbl calc_pressure(void)
Using this backend, calculate the pressure in Pa.
Definition: AbstractState.h:183
double smass(void)
Return the molar entropy in J/kg/K.
Definition: AbstractState.h:744
virtual CoolPropDbl calc_rhomolar(void)
Using this backend, get the molar density in mol/m^3.
Definition: AbstractState.h:394
CoolPropDbl d4alphar_dTau4(void)
Return the term .
Definition: AbstractState.h:1164
CoolPropDbl dalpha0_dDelta(void)
Return the term .
Definition: AbstractState.h:1048
virtual void set_reference_stateD(double T, double rhomolar, double hmolar0, double smolar0)
Set the reference state based on a thermodynamic state point specified by temperature and molar densi...
Definition: AbstractState.h:506
double Prandtl(void)
Return the Prandtl number (dimensionless)
Definition: AbstractState.h:1022
virtual CoolPropDbl calc_alpha0(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:236
double volumemass_excess(void)
Return the excess volume in m^3/kg.
Definition: AbstractState.h:790
virtual CoolPropDbl calc_rhomass_critical(void)
Using this backend, get the critical point mass density in kg/m^3 - Added for IF97Backend which is ma...
Definition: AbstractState.h:324
virtual CoolPropDbl calc_GWP20(void)
Using this backend, calculate the 20-year global warming potential (GWP)
Definition: AbstractState.h:266
virtual double get_binary_interaction_double(const std::string &CAS1, const std::string &CAS2, const std::string &parameter)
Get binary mixture double value (EXPERT USE ONLY!!!)
Definition: AbstractState.h:578
virtual CoolPropDbl calc_d2alphar_dDelta2(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:209
std::vector< double > M1
The determinant of the scaled matrix for the second criticality condition.
Definition: AbstractState.h:25
virtual std::vector< CoolPropDbl > calc_fugacity_coefficients()
Using this backend, calculate the fugacity in Pa.
Definition: AbstractState.h:189
virtual CoolPropDbl calc_cpmolar(void)
Using this backend, calculate the molar constant-pressure specific heat in J/mol/K.
Definition: AbstractState.h:153
virtual void set_binary_interaction_double(const std::size_t i, const std::size_t j, const std::string &parameter, const double value)
Set binary mixture floating point parameter (EXPERT USE ONLY!!!)
Definition: AbstractState.h:572
CoolPropDbl alpha0(void)
Return the term .
Definition: AbstractState.h:1043
virtual CoolPropDbl calc_gibbsmolar_residual(void)
Using this backend, calculate the residual molar Gibbs function in J/mol.
Definition: AbstractState.h:161
virtual std::string calc_name(void)
Using this backend, get the name of the fluid.
Definition: AbstractState.h:304
void unspecify_phase(void)
Unspecify the phase and go back to calculating it based on the inputs.
Definition: AbstractState.h:620
virtual CoolPropDbl calc_d4alphar_dDelta3_dTau(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:226
double p(void)
Return the pressure in Pa.
Definition: AbstractState.h:706
virtual CoolPropDbl calc_d3alphar_dDelta3(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:215
virtual CoolPropDbl calc_d4alphar_dDelta_dTau3(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:230
double umass(void)
Return the mass internal energy in J/kg.
Definition: AbstractState.h:752
virtual double get_binary_interaction_double(const std::size_t i, const std::size_t j, const std::string &parameter)
Get binary mixture double value (EXPERT USE ONLY!!!)
Definition: AbstractState.h:580
virtual CoolPropDbl calc_speed_sound(void)
Using this backend, calculate the speed of sound in m/s.
Definition: AbstractState.h:165
std::vector< double > mole_fractions_liquid_double(void)
Get the mole fractions of the equilibrium liquid phase (but as a double for use in SWIG wrapper) ...
Definition: AbstractState.h:524
CoolPropDbl d3alpha0_dDelta2_dTau(void)
Return the term .
Definition: AbstractState.h:1083
Twophase.
Definition: DataStructures.h:164
phases _phase
The key for the phase from CoolProp::phases enum.
Definition: AbstractState.h:75
virtual CoolPropDbl calc_T_critical(void)
Using this backend, get the critical point temperature in K.
Definition: AbstractState.h:314
CachedElement _hmolar_residual
Residual properties.
Definition: AbstractState.h:110
virtual CoolPropDbl calc_p_reducing(void)
Using this backend, get the reducing point pressure in Pa.
Definition: AbstractState.h:320
virtual const CoolProp::SimpleState & get_reducing_state()
Get the state that is used in the equation of state or mixture model to reduce the state...
Definition: AbstractState.h:601
double _rhomolar
Bulk values.
Definition: AbstractState.h:100
virtual CoolPropDbl calc_alphar(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:203
CoolPropDbl second_partial_deriv(parameters Of1, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2)
The second partial derivative in homogeneous phases.
Definition: AbstractState.h:864
std::vector< double > y
molar composition of the vapor phase
Definition: AbstractState.h:41
const CoolProp::PhaseEnvelopeData & get_phase_envelope_data()
After having calculated the phase envelope, return the phase envelope data.
Definition: AbstractState.h:988
virtual CoolPropDbl calc_d4alphar_dDelta4(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:224
This simple class holds the values for guesses for use in some solvers that have the ability to use g...
Definition: AbstractState.h:32
virtual CoolPropDbl calc_rhomolar_reducing(void)
Using this backend, get the reducing point molar density in mol/m^3.
Definition: AbstractState.h:326
const CoolProp::SimpleState & get_state(const std::string &state)
Get a desired state point - backend dependent.
Definition: AbstractState.h:604
virtual CoolPropDbl calc_dCvirial_dT(void)
Using this backend, calculate the derivative dC/dT.
Definition: AbstractState.h:299
CoolPropDbl alphar(void)
Return the term .
Definition: AbstractState.h:1094
CoolPropDbl d2alphar_dDelta2(void)
Return the term .
Definition: AbstractState.h:1109
virtual CoolPropDbl calc_p_triple(void)
Using this backend, get the triple point pressure in Pa.
Definition: AbstractState.h:311
double gibbsmass_excess(void)
Return the excess Gibbs energy in J/kg.
Definition: AbstractState.h:778
virtual CoolPropDbl calc_GWP100(void)
Using this backend, calculate the 100-year global warming potential (GWP)
Definition: AbstractState.h:268
Definition: AbstractState.h:1189
virtual void set_fluid_parameter_double(const size_t i, const std::string &parameter, const double value)
Set fluid parameter (currently the volume translation parameter for cubic)
Definition: AbstractState.h:588
virtual CoolPropDbl calc_fugacity(std::size_t i)
Using this backend, calculate the fugacity in Pa.
Definition: AbstractState.h:191
CachedElement _rhoLanc
Ancillary values.
Definition: AbstractState.h:116
virtual CoolPropDbl calc_p_critical(void)
Using this backend, get the critical point pressure in Pa.
Definition: AbstractState.h:318
CoolPropDbl d2alphar_dTau2(void)
Return the term .
Definition: AbstractState.h:1119
virtual CoolPropDbl calc_PIP(void)
Using this backend, calculate the phase identification parameter (PIP)
Definition: AbstractState.h:195
CoolPropDbl d3alphar_dDelta_dTau2(void)
Return the term .
Definition: AbstractState.h:1134
virtual CoolPropDbl calc_cvmolar(void)
Using this backend, calculate the molar constant-volume specific heat in J/mol/K. ...
Definition: AbstractState.h:157
virtual CoolPropDbl calc_d3alphar_dDelta2_dTau(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:217
virtual const CoolProp::SimpleState & calc_state(const std::string &state)
Using this backend, calculate a phase given by the state string.
Definition: AbstractState.h:367
virtual CoolPropDbl calc_fugacity_coefficient(std::size_t i)
Using this backend, calculate the fugacity coefficient (dimensionless)
Definition: AbstractState.h:187
virtual const std::vector< CoolPropDbl > get_mass_fractions(void)
Get the mass fractions of the fluid.
Definition: AbstractState.h:540
virtual CoolPropDbl calc_d4alphar_dTau4(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:232
Supercritical (p > pc, T > Tc)
Definition: DataStructures.h:159
This file contains flash routines in which the state is unknown, and a solver of some kind must be us...
Definition: AbstractState.h:19
CoolPropDbl d2alpha0_dDelta_dTau(void)
Return the term .
Definition: AbstractState.h:1063
virtual CoolPropDbl calc_fraction_max(void)
Get the maximum fraction (mole, mass, volume) for incompressible fluid.
Definition: AbstractState.h:378
virtual CoolPropDbl calc_smolar_residual(void)
Using this backend, calculate the residual molar entropy in J/mol/K.
Definition: AbstractState.h:149
void viscosity_contributions(CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical)
Return the viscosity contributions, each in Pa-s.
Definition: AbstractState.h:1014
CoolPropDbl dalphar_dTau(void)
Return the term .
Definition: AbstractState.h:1104
virtual void set_reference_stateS(const std::string &reference_state)
Set the reference state based on a string representation.
Definition: AbstractState.h:497
virtual CoolPropDbl calc_d3alpha0_dDelta2_dTau(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:250
CoolPropDbl d4alphar_dDelta4(void)
Return the term .
Definition: AbstractState.h:1144
parameters
Define some constants that will be used throughout These are constants for the input and output para...
Definition: DataStructures.h:49
virtual CoolPropDbl calc_gibbsmolar(void)
Using this backend, calculate the molar Gibbs function in J/mol.
Definition: AbstractState.h:159
virtual CoolPropDbl calc_d3alpha0_dDelta_dTau2(void)
Using this backend, calculate the ideal-gas Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:252
virtual phases calc_phase(void)
Using this backend, calculate the phase.
Definition: AbstractState.h:358
void register_backend(const backend_families &bf, shared_ptr< AbstractStateGenerator > gen)
Register a backend in the backend library (statically defined in AbstractState.cpp and not publicly a...
Definition: AbstractState.cpp:52
Definition: CachedElement.h:32
double first_two_phase_deriv(parameters Of, parameters Wrt, parameters Constant)
Calculate the first "two-phase" derivative as described by Thorade and Sadaat, EAS, 2013.
Definition: AbstractState.h:927
virtual CoolPropDbl calc_dalphar_dTau(void)
Using this backend, calculate the residual Helmholtz energy term (dimensionless) ...
Definition: AbstractState.h:207
virtual CoolPropDbl calc_smolar(void)
Using this backend, calculate the molar entropy in J/mol/K.
Definition: AbstractState.h:147
virtual CoolPropDbl calc_rhomolar_critical(void)
Using this backend, get the critical point molar density in mol/m^3.
Definition: AbstractState.h:322
virtual CoolPropDbl calc_isothermal_compressibility(void)
Using this backend, calculate the isothermal compressibility in 1/Pa.
Definition: AbstractState.h:167
double second_two_phase_deriv(parameters Of, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2)
Calculate the second "two-phase" derivative as described by Thorade and Sadaat, EAS, 2013.
Definition: AbstractState.h:947
CoolPropDbl second_saturation_deriv(parameters Of1, parameters Wrt1, parameters Wrt2)
The second partial derivative along the saturation curve.
Definition: AbstractState.h:906
virtual CoolPropDbl calc_surface_tension(void)
Using this backend, calculate the surface tension in N/m.
Definition: AbstractState.h:177
virtual CoolPropDbl calc_dipole_moment(void)
Using this backend, calculate the dipole moment in C-m (1 D = 3.33564e-30 C-m)
Definition: AbstractState.h:280
void ideal_curve(const std::string &type, std::vector< double > &T, std::vector< double > &p)
Calculate an ideal curve for a pure fluid.
Definition: AbstractState.h:829
double dipole_moment()
Return the dipole moment in C-m (1 D = 3.33564e-30 C-m)
Definition: AbstractState.h:684
Definition: DataStructures.h:16
SpinodalData get_spinodal_data()
Get the data from the spinodal curve constructed in the call to build_spinodal()
Definition: AbstractState.h:638
std::vector< CoolPropDbl > mole_fractions_vapor(void)
Get the mole fractions of the equilibrium vapor phase.
Definition: AbstractState.h:530