CoolProp  6.6.0
An open-source fluid property and humid air property database
AbstractState.h
Go to the documentation of this file.
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"
16 
17 #include <numeric>
18 
19 namespace CoolProp {
20 
24 {
25  public:
26  std::vector<double> tau,
28  M1;
29 };
30 
34 {
35  public:
36  double T,
37  p,
43  std::vector<double> x,
44  y;
46  clear();
47  };
48  void clear() {
49  T = _HUGE;
50  p = _HUGE;
51  rhomolar = _HUGE;
52  hmolar = _HUGE;
53  smolar = _HUGE;
54  rhomolar_liq = _HUGE;
55  rhomolar_vap = _HUGE;
56  x.clear(), y.clear();
57  }
58 };
59 
61 
79 {
80  protected:
85 
86  bool isSupercriticalPhase(void) {
87  return (this->_phase == iphase_supercritical || this->_phase == iphase_supercritical_liquid || this->_phase == iphase_supercritical_gas);
88  }
89 
90  bool isHomogeneousPhase(void) {
91  return (this->_phase == iphase_liquid || this->_phase == iphase_gas || isSupercriticalPhase() || this->_phase == iphase_critical_point);
92  }
93 
94  bool isTwoPhase(void) {
95  return (this->_phase == iphase_twophase);
96  }
97 
100 
103 
106 
108  double _rhomolar, _T, _p, _Q, _R;
109 
111 
114 
116 
119 
122 
125 
127 
130 
136 
138 
141 
142  // ----------------------------------------
143  // Property accessors to be optionally implemented by the backend
144  // for properties that are not always calculated
145  // ----------------------------------------
147  virtual CoolPropDbl calc_hmolar(void) {
148  throw NotImplementedError("calc_hmolar is not implemented for this backend");
149  };
152  throw NotImplementedError("calc_hmolar_residual is not implemented for this backend");
153  };
155  virtual CoolPropDbl calc_smolar(void) {
156  throw NotImplementedError("calc_smolar is not implemented for this backend");
157  };
160  throw NotImplementedError("calc_smolar_residual is not implemented for this backend");
161  };
163  virtual CoolPropDbl calc_umolar(void) {
164  throw NotImplementedError("calc_umolar is not implemented for this backend");
165  };
167  virtual CoolPropDbl calc_cpmolar(void) {
168  throw NotImplementedError("calc_cpmolar is not implemented for this backend");
169  };
172  throw NotImplementedError("calc_cpmolar_idealgas is not implemented for this backend");
173  };
175  virtual CoolPropDbl calc_cvmolar(void) {
176  throw NotImplementedError("calc_cvmolar is not implemented for this backend");
177  };
179  virtual CoolPropDbl calc_gibbsmolar(void) {
180  throw NotImplementedError("calc_gibbsmolar is not implemented for this backend");
181  };
184  throw NotImplementedError("calc_gibbsmolar_residual is not implemented for this backend");
185  };
188  throw NotImplementedError("calc_helmholtzmolar is not implemented for this backend");
189  };
192  throw NotImplementedError("calc_speed_sound is not implemented for this backend");
193  };
196  throw NotImplementedError("calc_isothermal_compressibility is not implemented for this backend");
197  };
200  throw NotImplementedError("calc_isobaric_expansion_coefficient is not implemented for this backend");
201  };
204  throw NotImplementedError("calc_isentropic_expansion_coefficient is not implemented for this backend");
205  };
207  virtual CoolPropDbl calc_viscosity(void) {
208  throw NotImplementedError("calc_viscosity is not implemented for this backend");
209  };
212  throw NotImplementedError("calc_conductivity is not implemented for this backend");
213  };
216  throw NotImplementedError("calc_surface_tension is not implemented for this backend");
217  };
219  virtual CoolPropDbl calc_molar_mass(void) {
220  throw NotImplementedError("calc_molar_mass is not implemented for this backend");
221  };
224  throw NotImplementedError("calc_acentric_factor is not implemented for this backend");
225  };
227  virtual CoolPropDbl calc_pressure(void) {
228  throw NotImplementedError("calc_pressure is not implemented for this backend");
229  };
232  throw NotImplementedError("calc_gas_constant is not implemented for this backend");
233  };
235  virtual CoolPropDbl calc_fugacity_coefficient(std::size_t i) {
236  throw NotImplementedError("calc_fugacity_coefficient is not implemented for this backend");
237  };
239  virtual std::vector<CoolPropDbl> calc_fugacity_coefficients() {
240  throw NotImplementedError("calc_fugacity_coefficients is not implemented for this backend");
241  };
243  virtual CoolPropDbl calc_fugacity(std::size_t i) {
244  throw NotImplementedError("calc_fugacity is not implemented for this backend");
245  };
247  virtual CoolPropDbl calc_chemical_potential(std::size_t i) {
248  throw NotImplementedError("calc_chemical_potential is not implemented for this backend");
249  };
251  virtual CoolPropDbl calc_PIP(void) {
252  throw NotImplementedError("calc_PIP is not implemented for this backend");
253  };
254 
255  // Excess properties
257  virtual void calc_excess_properties(void) {
258  throw NotImplementedError("calc_excess_properties is not implemented for this backend");
259  };
260 
261  // Derivatives of residual helmholtz energy
263  virtual CoolPropDbl calc_alphar(void) {
264  throw NotImplementedError("calc_alphar is not implemented for this backend");
265  };
268  throw NotImplementedError("calc_dalphar_dDelta is not implemented for this backend");
269  };
272  throw NotImplementedError("calc_dalphar_dTau is not implemented for this backend");
273  };
276  throw NotImplementedError("calc_d2alphar_dDelta2 is not implemented for this backend");
277  };
280  throw NotImplementedError("calc_d2alphar_dDelta_dTau is not implemented for this backend");
281  };
284  throw NotImplementedError("calc_d2alphar_dTau2 is not implemented for this backend");
285  };
288  throw NotImplementedError("calc_d3alphar_dDelta3 is not implemented for this backend");
289  };
292  throw NotImplementedError("calc_d3alphar_dDelta2_dTau is not implemented for this backend");
293  };
296  throw NotImplementedError("calc_d3alphar_dDelta_dTau2 is not implemented for this backend");
297  };
300  throw NotImplementedError("calc_d3alphar_dTau3 is not implemented for this backend");
301  };
302 
305  throw NotImplementedError("calc_d4alphar_dDelta4 is not implemented for this backend");
306  };
309  throw NotImplementedError("calc_d4alphar_dDelta3_dTau is not implemented for this backend");
310  };
313  throw NotImplementedError("calc_d4alphar_dDelta2_dTau2 is not implemented for this backend");
314  };
317  throw NotImplementedError("calc_d4alphar_dDelta_dTau3 is not implemented for this backend");
318  };
321  throw NotImplementedError("calc_d4alphar_dTau4 is not implemented for this backend");
322  };
323 
324  // Derivatives of ideal-gas helmholtz energy
326  virtual CoolPropDbl calc_alpha0(void) {
327  throw NotImplementedError("calc_alpha0 is not implemented for this backend");
328  };
331  throw NotImplementedError("calc_dalpha0_dDelta is not implemented for this backend");
332  };
335  throw NotImplementedError("calc_dalpha0_dTau is not implemented for this backend");
336  };
339  throw NotImplementedError("calc_d2alpha0_dDelta_dTau is not implemented for this backend");
340  };
343  throw NotImplementedError("calc_d2alpha0_dDelta2 is not implemented for this backend");
344  };
347  throw NotImplementedError("calc_d2alpha0_dTau2 is not implemented for this backend");
348  };
351  throw NotImplementedError("calc_d3alpha0_dDelta3 is not implemented for this backend");
352  };
355  throw NotImplementedError("calc_d3alpha0_dDelta2_dTau is not implemented for this backend");
356  };
359  throw NotImplementedError("calc_d3alpha0_dDelta_dTau2 is not implemented for this backend");
360  };
363  throw NotImplementedError("calc_d3alpha0_dTau3 is not implemented for this backend");
364  };
365 
366  virtual void calc_reducing_state(void) {
367  throw NotImplementedError("calc_reducing_state is not implemented for this backend");
368  };
369 
371  virtual CoolPropDbl calc_Tmax(void) {
372  throw NotImplementedError("calc_Tmax is not implemented for this backend");
373  };
375  virtual CoolPropDbl calc_Tmin(void) {
376  throw NotImplementedError("calc_Tmin is not implemented for this backend");
377  };
379  virtual CoolPropDbl calc_pmax(void) {
380  throw NotImplementedError("calc_pmax is not implemented for this backend");
381  };
382 
384  virtual CoolPropDbl calc_GWP20(void) {
385  throw NotImplementedError("calc_GWP20 is not implemented for this backend");
386  };
388  virtual CoolPropDbl calc_GWP100(void) {
389  throw NotImplementedError("calc_GWP100 is not implemented for this backend");
390  };
392  virtual CoolPropDbl calc_GWP500(void) {
393  throw NotImplementedError("calc_GWP500 is not implemented for this backend");
394  };
396  virtual CoolPropDbl calc_ODP(void) {
397  throw NotImplementedError("calc_ODP is not implemented for this backend");
398  };
401  throw NotImplementedError("calc_flame_hazard is not implemented for this backend");
402  };
405  throw NotImplementedError("calc_health_hazard is not implemented for this backend");
406  };
409  throw NotImplementedError("calc_physical_hazard is not implemented for this backend");
410  };
413  throw NotImplementedError("calc_dipole_moment is not implemented for this backend");
414  };
415 
419  virtual CoolPropDbl calc_second_partial_deriv(parameters Of1, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2);
420 
423  throw NotImplementedError("calc_reduced_density is not implemented for this backend");
424  };
427  throw NotImplementedError("calc_reciprocal_reduced_temperature is not implemented for this backend");
428  };
429 
431  virtual CoolPropDbl calc_Bvirial(void) {
432  throw NotImplementedError("calc_Bvirial is not implemented for this backend");
433  };
435  virtual CoolPropDbl calc_Cvirial(void) {
436  throw NotImplementedError("calc_Cvirial is not implemented for this backend");
437  };
440  throw NotImplementedError("calc_dBvirial_dT is not implemented for this backend");
441  };
444  throw NotImplementedError("calc_dCvirial_dT is not implemented for this backend");
445  };
448  throw NotImplementedError("calc_compressibility_factor is not implemented for this backend");
449  };
450 
452  virtual std::string calc_name(void) {
453  throw NotImplementedError("calc_name is not implemented for this backend");
454  };
456  virtual std::string calc_description(void) {
457  throw NotImplementedError("calc_description is not implemented for this backend");
458  };
459 
461  virtual CoolPropDbl calc_Ttriple(void) {
462  throw NotImplementedError("calc_Ttriple is not implemented for this backend");
463  };
465  virtual CoolPropDbl calc_p_triple(void) {
466  throw NotImplementedError("calc_p_triple is not implemented for this backend");
467  };
468 
470  virtual CoolPropDbl calc_T_critical(void) {
471  throw NotImplementedError("calc_T_critical is not implemented for this backend");
472  };
474  virtual CoolPropDbl calc_T_reducing(void) {
475  throw NotImplementedError("calc_T_reducing is not implemented for this backend");
476  };
478  virtual CoolPropDbl calc_p_critical(void) {
479  throw NotImplementedError("calc_p_critical is not implemented for this backend");
480  };
482  virtual CoolPropDbl calc_p_reducing(void) {
483  throw NotImplementedError("calc_p_reducing is not implemented for this backend");
484  };
487  throw NotImplementedError("calc_rhomolar_critical is not implemented for this backend");
488  };
491  throw NotImplementedError("calc_rhomass_critical is not implemented for this backend");
492  };
495  throw NotImplementedError("calc_rhomolar_reducing is not implemented for this backend");
496  };
498  virtual void calc_phase_envelope(const std::string& type) {
499  throw NotImplementedError("calc_phase_envelope is not implemented for this backend");
500  };
502  virtual CoolPropDbl calc_rhomass(void) {
503  return rhomolar() * molar_mass();
504  }
505  virtual CoolPropDbl calc_hmass(void) {
506  return hmolar() / molar_mass();
507  }
509  return hmolar_excess() / molar_mass();
510  }
511  virtual CoolPropDbl calc_smass(void) {
512  return smolar() / molar_mass();
513  }
515  return smolar_excess() / molar_mass();
516  }
517  virtual CoolPropDbl calc_cpmass(void) {
518  return cpmolar() / molar_mass();
519  }
520  virtual CoolPropDbl calc_cp0mass(void) {
521  return cp0molar() / molar_mass();
522  }
523  virtual CoolPropDbl calc_cvmass(void) {
524  return cvmolar() / molar_mass();
525  }
526  virtual CoolPropDbl calc_umass(void) {
527  return umolar() / molar_mass();
528  }
530  return umolar_excess() / molar_mass();
531  }
532  virtual CoolPropDbl calc_gibbsmass(void) {
533  return gibbsmolar() / molar_mass();
534  }
536  return gibbsmolar_excess() / molar_mass();
537  }
539  return helmholtzmolar() / molar_mass();
540  }
542  return helmholtzmolar_excess() / molar_mass();
543  }
545  return volumemolar_excess() / molar_mass();
546  }
547 
549  virtual void update_states(void) {
550  throw NotImplementedError("This backend does not implement update_states function");
551  };
552 
553  virtual CoolPropDbl calc_melting_line(int param, int given, CoolPropDbl value) {
554  throw NotImplementedError("This backend does not implement calc_melting_line function");
555  };
556 
561  virtual CoolPropDbl calc_saturation_ancillary(parameters param, int Q, parameters given, double value) {
562  throw NotImplementedError("This backend does not implement calc_saturation_ancillary");
563  };
564 
566  virtual phases calc_phase(void) {
567  throw NotImplementedError("This backend does not implement calc_phase function");
568  };
571  throw NotImplementedError("This backend does not implement calc_specify_phase function");
572  };
574  virtual void calc_unspecify_phase(void) {
575  throw NotImplementedError("This backend does not implement calc_unspecify_phase function");
576  };
578  virtual std::vector<std::string> calc_fluid_names(void) {
579  throw NotImplementedError("This backend does not implement calc_fluid_names function");
580  };
583  virtual const CoolProp::SimpleState& calc_state(const std::string& state) {
584  throw NotImplementedError("calc_state is not implemented for this backend");
585  };
586 
588  throw NotImplementedError("calc_phase_envelope_data is not implemented for this backend");
589  };
590 
591  virtual std::vector<CoolPropDbl> calc_mole_fractions_liquid(void) {
592  throw NotImplementedError("calc_mole_fractions_liquid is not implemented for this backend");
593  };
594  virtual std::vector<CoolPropDbl> calc_mole_fractions_vapor(void) {
595  throw NotImplementedError("calc_mole_fractions_vapor is not implemented for this backend");
596  };
597  virtual const std::vector<CoolPropDbl> calc_mass_fractions(void) {
598  throw NotImplementedError("calc_mass_fractions is not implemented for this backend");
599  };
600 
603  throw NotImplementedError("calc_fraction_min is not implemented for this backend");
604  };
607  throw NotImplementedError("calc_fraction_max is not implemented for this backend");
608  };
609  virtual CoolPropDbl calc_T_freeze(void) {
610  throw NotImplementedError("calc_T_freeze is not implemented for this backend");
611  };
612 
614  throw NotImplementedError("calc_first_saturation_deriv is not implemented for this backend");
615  };
617  throw NotImplementedError("calc_second_saturation_deriv is not implemented for this backend");
618  };
620  throw NotImplementedError("calc_first_two_phase_deriv is not implemented for this backend");
621  };
623  throw NotImplementedError("calc_second_two_phase_deriv is not implemented for this backend");
624  };
626  throw NotImplementedError("calc_first_two_phase_deriv_splined is not implemented for this backend");
627  };
628 
630  throw NotImplementedError("calc_saturated_liquid_keyed_output is not implemented for this backend");
631  };
633  throw NotImplementedError("calc_saturated_vapor_keyed_output is not implemented for this backend");
634  };
635  virtual void calc_ideal_curve(const std::string& type, std::vector<double>& T, std::vector<double>& p) {
636  throw NotImplementedError("calc_ideal_curve is not implemented for this backend");
637  };
638 
640  virtual CoolPropDbl calc_T(void) {
641  return _T;
642  }
644  virtual CoolPropDbl calc_rhomolar(void) {
645  return _rhomolar;
646  }
647 
649  virtual double calc_tangent_plane_distance(const double T, const double p, const std::vector<double>& w, const double rhomolar_guess) {
650  throw NotImplementedError("calc_tangent_plane_distance is not implemented for this backend");
651  };
652 
654  virtual void calc_true_critical_point(double& T, double& rho) {
655  throw NotImplementedError("calc_true_critical_point is not implemented for this backend");
656  };
657 
658  virtual void calc_conformal_state(const std::string& reference_fluid, CoolPropDbl& T, CoolPropDbl& rhomolar) {
659  throw NotImplementedError("calc_conformal_state is not implemented for this backend");
660  };
661 
662  virtual void calc_viscosity_contributions(CoolPropDbl& dilute, CoolPropDbl& initial_density, CoolPropDbl& residual, CoolPropDbl& critical) {
663  throw NotImplementedError("calc_viscosity_contributions is not implemented for this backend");
664  };
665  virtual void calc_conductivity_contributions(CoolPropDbl& dilute, CoolPropDbl& initial_density, CoolPropDbl& residual, CoolPropDbl& critical) {
666  throw NotImplementedError("calc_conductivity_contributions is not implemented for this backend");
667  };
668  virtual std::vector<CriticalState> calc_all_critical_points(void) {
669  throw NotImplementedError("calc_all_critical_points is not implemented for this backend");
670  };
671  virtual void calc_build_spinodal() {
672  throw NotImplementedError("calc_build_spinodal is not implemented for this backend");
673  };
675  throw NotImplementedError("calc_get_spinodal_data is not implemented for this backend");
676  };
677  virtual void calc_criticality_contour_values(double& L1star, double& M1star) {
678  throw NotImplementedError("calc_criticality_contour_values is not implemented for this backend");
679  };
680 
682  virtual void mass_to_molar_inputs(CoolProp::input_pairs& input_pair, CoolPropDbl& value1, CoolPropDbl& value2);
683 
685  virtual void calc_change_EOS(const std::size_t i, const std::string& EOS_name) {
686  throw NotImplementedError("calc_change_EOS is not implemented for this backend");
687  };
688 
689  public:
691  clear();
692  }
693  virtual ~AbstractState(){};
694 
696 
702  static AbstractState* factory(const std::string& backend, const std::string& fluid_names) {
703  return factory(backend, strsplit(fluid_names, '&'));
704  };
705 
729  static AbstractState* factory(const std::string& backend, const std::vector<std::string>& fluid_names);
730 
733  _T = T;
734  }
735 
740  virtual std::string backend_name(void) = 0;
741 
742  // The derived classes must implement this function to define whether they use mole fractions (true) or mass fractions (false)
743  virtual bool using_mole_fractions(void) = 0;
744  virtual bool using_mass_fractions(void) = 0;
745  virtual bool using_volu_fractions(void) = 0;
746 
747  virtual void set_mole_fractions(const std::vector<CoolPropDbl>& mole_fractions) = 0;
748  virtual void set_mass_fractions(const std::vector<CoolPropDbl>& mass_fractions) = 0;
749  virtual void set_volu_fractions(const std::vector<CoolPropDbl>& mass_fractions) {
750  throw NotImplementedError("Volume composition has not been implemented.");
751  }
752 
772  virtual void set_reference_stateS(const std::string& reference_state) {
773  throw NotImplementedError(
774  "Setting reference state has not been implemented for this backend. Try using CoolProp::set_reference_stateD instead.");
775  }
776 
782  virtual void set_reference_stateD(double T, double rhomolar, double hmolar0, double smolar0) {
783  throw NotImplementedError(
784  "Setting reference state has not been implemented for this backend. Try using CoolProp::set_reference_stateD instead.");
785  }
786 
787 #ifndef COOLPROPDBL_MAPS_TO_DOUBLE
788  void set_mole_fractions(const std::vector<double>& mole_fractions) {
789  set_mole_fractions(std::vector<CoolPropDbl>(mole_fractions.begin(), mole_fractions.end()));
790  };
791  void set_mass_fractions(const std::vector<double>& mass_fractions) {
792  set_mass_fractions(std::vector<CoolPropDbl>(mass_fractions.begin(), mass_fractions.end()));
793  };
794  void set_volu_fractions(const std::vector<double>& volu_fractions) {
795  set_volu_fractions(std::vector<CoolPropDbl>(volu_fractions.begin(), volu_fractions.end()));
796  };
797 #endif
798 
799 #ifdef EMSCRIPTEN
800  void set_mole_fractions_double(const std::vector<double>& mole_fractions) {
801  set_mole_fractions(std::vector<CoolPropDbl>(mole_fractions.begin(), mole_fractions.end()));
802  };
803 #endif
804 
806  std::vector<CoolPropDbl> mole_fractions_liquid(void) {
808  };
810  std::vector<double> mole_fractions_liquid_double(void) {
811  std::vector<CoolPropDbl> x = calc_mole_fractions_liquid();
812  return std::vector<double>(x.begin(), x.end());
813  };
814 
816  std::vector<CoolPropDbl> mole_fractions_vapor(void) {
817  return calc_mole_fractions_vapor();
818  };
820  std::vector<double> mole_fractions_vapor_double(void) {
821  std::vector<CoolPropDbl> y = calc_mole_fractions_vapor();
822  return std::vector<double>(y.begin(), y.end());
823  };
824 
826  virtual const std::vector<CoolPropDbl>& get_mole_fractions(void) = 0;
828  virtual const std::vector<CoolPropDbl> get_mass_fractions(void) {
829  return this->calc_mass_fractions();
830  };
831 
833  virtual void update(CoolProp::input_pairs input_pair, double Value1, double Value2) = 0;
834 
837  virtual void update_with_guesses(CoolProp::input_pairs input_pair, double Value1, double Value2, const GuessesStructure& guesses) {
838  throw NotImplementedError("update_with_guesses is not implemented for this backend");
839  };
840 
844  virtual bool available_in_high_level(void) {
845  return true;
846  }
847 
849  virtual std::string fluid_param_string(const std::string&) {
850  throw NotImplementedError("fluid_param_string has not been implemented for this backend");
851  }
852 
854  std::vector<std::string> fluid_names(void);
855 
860  virtual const double get_fluid_constant(std::size_t i, parameters param) const {
861  throw NotImplementedError("get_fluid_constant is not implemented for this backend");
862  };
863  ;
864 
866  virtual void set_binary_interaction_double(const std::string& CAS1, const std::string& CAS2, const std::string& parameter, const double value) {
867  throw NotImplementedError("set_binary_interaction_double is not implemented for this backend");
868  };
870  virtual void set_binary_interaction_double(const std::size_t i, const std::size_t j, const std::string& parameter, const double value) {
871  throw NotImplementedError("set_binary_interaction_double is not implemented for this backend");
872  };
874  virtual void set_binary_interaction_string(const std::string& CAS1, const std::string& CAS2, const std::string& parameter,
875  const std::string& value) {
876  throw NotImplementedError("set_binary_interaction_string is not implemented for this backend");
877  };
879  virtual void set_binary_interaction_string(const std::size_t i, const std::size_t j, const std::string& parameter, const std::string& value) {
880  throw NotImplementedError("set_binary_interaction_string is not implemented for this backend");
881  };
883  virtual double get_binary_interaction_double(const std::string& CAS1, const std::string& CAS2, const std::string& parameter) {
884  throw NotImplementedError("get_binary_interaction_double is not implemented for this backend");
885  };
887  virtual double get_binary_interaction_double(const std::size_t i, const std::size_t j, const std::string& parameter) {
888  throw NotImplementedError("get_binary_interaction_double is not implemented for this backend");
889  };
891  virtual std::string get_binary_interaction_string(const std::string& CAS1, const std::string& CAS2, const std::string& parameter) {
892  throw NotImplementedError("get_binary_interaction_string is not implemented for this backend");
893  };
895  virtual void apply_simple_mixing_rule(std::size_t i, std::size_t j, const std::string& model) {
896  throw NotImplementedError("apply_simple_mixing_rule is not implemented for this backend");
897  };
899  virtual void set_cubic_alpha_C(const size_t i, const std::string& parameter, const double c1, const double c2, const double c3) {
900  throw ValueError("set_cubic_alpha_C only defined for cubic backends");
901  };
903  virtual void set_fluid_parameter_double(const size_t i, const std::string& parameter, const double value) {
904  throw ValueError("set_fluid_parameter_double only defined for cubic backends");
905  };
907  virtual double get_fluid_parameter_double(const size_t i, const std::string& parameter) {
908  throw ValueError("get_fluid_parameter_double only defined for cubic backends");
909  };
910 
912  virtual bool clear();
914  virtual bool clear_comp_change();
915 
920  return _reducing;
921  };
922 
924  const CoolProp::SimpleState& get_state(const std::string& state) {
925  return calc_state(state);
926  };
927 
929  double Tmin(void);
931  double Tmax(void);
933  double pmax(void);
935  double Ttriple(void);
936 
938  phases phase(void) {
939  return calc_phase();
940  };
944  };
946  void unspecify_phase(void) {
948  };
949 
951  double T_critical(void);
953  double p_critical(void);
955  double rhomolar_critical(void);
957  double rhomass_critical(void);
958 
960  std::vector<CriticalState> all_critical_points(void) {
961  return calc_all_critical_points();
962  };
963 
965  void build_spinodal() {
967  };
968 
971  return calc_get_spinodal_data();
972  };
973 
975  void criticality_contour_values(double& L1star, double& M1star) {
976  return calc_criticality_contour_values(L1star, M1star);
977  }
978 
1002  double tangent_plane_distance(const double T, const double p, const std::vector<double>& w, const double rhomolar_guess = -1) {
1003  return calc_tangent_plane_distance(T, p, w, rhomolar_guess);
1004  };
1005 
1007  double T_reducing(void);
1009  double rhomolar_reducing(void);
1011  double rhomass_reducing(void);
1012 
1014  double p_triple(void);
1015 
1017  std::string name() {
1018  return calc_name();
1019  };
1021  std::string description() {
1022  return calc_description();
1023  };
1024 
1026  double dipole_moment() {
1027  return calc_dipole_moment();
1028  }
1029 
1030  // ----------------------------------------
1031  // Bulk properties - temperature and density are directly calculated every time
1032  // All other parameters are calculated on an as-needed basis
1033  // ----------------------------------------
1035  double keyed_output(parameters key);
1037  double trivial_keyed_output(parameters key);
1041  };
1045  };
1046 
1048  double T(void) {
1049  return calc_T();
1050  };
1052  double rhomolar(void) {
1053  return calc_rhomolar();
1054  };
1056  double rhomass(void) {
1057  return calc_rhomass();
1058  };
1060  double p(void) {
1061  return _p;
1062  };
1064  double Q(void) {
1065  return _Q;
1066  };
1068  double tau(void);
1070  double delta(void);
1072  double molar_mass(void);
1074  double acentric_factor(void);
1076  double gas_constant(void);
1078  double Bvirial(void);
1080  double dBvirial_dT(void);
1082  double Cvirial(void);
1084  double dCvirial_dT(void);
1086  double compressibility_factor(void);
1088  double hmolar(void);
1090  double hmolar_residual(void);
1092  double hmass(void) {
1093  return calc_hmass();
1094  };
1096  double hmolar_excess(void);
1098  double hmass_excess(void) {
1099  return calc_hmass_excess();
1100  };
1102  double smolar(void);
1104  double smolar_residual(void);
1106  double smass(void) {
1107  return calc_smass();
1108  };
1110  double smolar_excess(void);
1112  double smass_excess(void) {
1113  return calc_smass_excess();
1114  };
1116  double umolar(void);
1118  double umass(void) {
1119  return calc_umass();
1120  };
1122  double umolar_excess(void);
1124  double umass_excess(void) {
1125  return calc_umass_excess();
1126  };
1128  double cpmolar(void);
1130  double cpmass(void) {
1131  return calc_cpmass();
1132  };
1134  double cp0molar(void);
1136  double cp0mass(void) {
1137  return calc_cp0mass();
1138  };
1140  double cvmolar(void);
1142  double cvmass(void) {
1143  return calc_cvmass();
1144  };
1146  double gibbsmolar(void);
1148  double gibbsmolar_residual(void);
1150  double gibbsmass(void) {
1151  return calc_gibbsmass();
1152  };
1154  double gibbsmolar_excess(void);
1156  double gibbsmass_excess(void) {
1157  return calc_gibbsmass_excess();
1158  };
1160  double helmholtzmolar(void);
1162  double helmholtzmass(void) {
1163  return calc_helmholtzmass();
1164  };
1166  double helmholtzmolar_excess(void);
1168  double helmholtzmass_excess(void) {
1169  return calc_helmholtzmass_excess();
1170  };
1172  double volumemolar_excess(void);
1174  double volumemass_excess(void) {
1175  return calc_volumemass_excess();
1176  };
1178  double speed_sound(void);
1180  double isothermal_compressibility(void);
1182  double isobaric_expansion_coefficient(void);
1184  double isentropic_expansion_coefficient(void);
1186  double fugacity_coefficient(std::size_t i);
1188  std::vector<double> fugacity_coefficients();
1190  double fugacity(std::size_t i);
1192  double chemical_potential(std::size_t i);
1203  double PIP() {
1204  return calc_PIP();
1205  };
1206 
1208  void true_critical_point(double& T, double& rho) {
1210  }
1211 
1219  void ideal_curve(const std::string& type, std::vector<double>& T, std::vector<double>& p) {
1220  calc_ideal_curve(type, T, p);
1221  };
1222 
1223  // ----------------------------------------
1224  // Partial derivatives
1225  // ----------------------------------------
1226 
1232  return calc_first_partial_deriv(Of, Wrt, Constant);
1233  };
1234 
1259  return calc_second_partial_deriv(Of1, Wrt1, Constant1, Wrt2, Constant2);
1260  };
1261 
1284  return calc_first_saturation_deriv(Of1, Wrt1);
1285  };
1286 
1305  return calc_second_saturation_deriv(Of1, Wrt1, Wrt2);
1306  };
1307 
1328  return calc_first_two_phase_deriv(Of, Wrt, Constant);
1329  };
1330 
1347  double second_two_phase_deriv(parameters Of, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2) {
1348  return calc_second_two_phase_deriv(Of, Wrt1, Constant1, Wrt2, Constant2);
1349  };
1350 
1371  double first_two_phase_deriv_splined(parameters Of, parameters Wrt, parameters Constant, double x_end) {
1372  return calc_first_two_phase_deriv_splined(Of, Wrt, Constant, x_end);
1373  };
1374 
1375  // ----------------------------------------
1376  // Phase envelope for mixtures
1377  // ----------------------------------------
1378 
1384  void build_phase_envelope(const std::string& type = "");
1389  return calc_phase_envelope_data();
1390  };
1391 
1392  // ----------------------------------------
1393  // Ancillary equations
1394  // ----------------------------------------
1395 
1397  virtual bool has_melting_line(void) {
1398  return false;
1399  };
1404  double melting_line(int param, int given, double value);
1410  double saturation_ancillary(parameters param, int Q, parameters given, double value);
1411 
1412  // ----------------------------------------
1413  // Transport properties
1414  // ----------------------------------------
1416  double viscosity(void);
1418  void viscosity_contributions(CoolPropDbl& dilute, CoolPropDbl& initial_density, CoolPropDbl& residual, CoolPropDbl& critical) {
1419  calc_viscosity_contributions(dilute, initial_density, residual, critical);
1420  };
1422  double conductivity(void);
1424  void conductivity_contributions(CoolPropDbl& dilute, CoolPropDbl& initial_density, CoolPropDbl& residual, CoolPropDbl& critical) {
1425  calc_conductivity_contributions(dilute, initial_density, residual, critical);
1426  };
1428  double surface_tension(void);
1430  double Prandtl(void) {
1431  return cpmass() * viscosity() / conductivity();
1432  };
1439  void conformal_state(const std::string& reference_fluid, CoolPropDbl& T, CoolPropDbl& rhomolar) {
1440  return calc_conformal_state(reference_fluid, T, rhomolar);
1441  };
1442 
1447  void change_EOS(const std::size_t i, const std::string& EOS_name) {
1448  calc_change_EOS(i, EOS_name);
1449  }
1450 
1451  // ----------------------------------------
1452  // Helmholtz energy and derivatives
1453  // ----------------------------------------
1456  if (!_alpha0) _alpha0 = calc_alpha0();
1457  return _alpha0;
1458  };
1462  return _dalpha0_dDelta;
1463  };
1467  return _dalpha0_dTau;
1468  };
1472  return _d2alpha0_dDelta2;
1473  };
1477  return _d2alpha0_dDelta_dTau;
1478  };
1482  return _d2alpha0_dTau2;
1483  };
1487  return _d3alpha0_dTau3;
1488  };
1492  return _d3alpha0_dDelta_dTau2;
1493  };
1497  return _d3alpha0_dDelta2_dTau;
1498  };
1502  return _d3alpha0_dDelta3;
1503  };
1504 
1507  if (!_alphar) _alphar = calc_alphar();
1508  return _alphar;
1509  };
1513  return _dalphar_dDelta;
1514  };
1518  return _dalphar_dTau;
1519  };
1523  return _d2alphar_dDelta2;
1524  };
1528  return _d2alphar_dDelta_dTau;
1529  };
1533  return _d2alphar_dTau2;
1534  };
1538  return _d3alphar_dDelta3;
1539  };
1543  return _d3alphar_dDelta2_dTau;
1544  };
1548  return _d3alphar_dDelta_dTau2;
1549  };
1553  return _d3alphar_dTau3;
1554  };
1558  return _d4alphar_dDelta4;
1559  };
1563  return _d4alphar_dDelta3_dTau;
1564  };
1568  return _d4alphar_dDelta2_dTau2;
1569  };
1573  return _d4alphar_dDelta_dTau3;
1574  };
1578  return _d4alphar_dTau4;
1579  };
1580 };
1581 
1590 {
1591  public:
1592  virtual AbstractState* get_AbstractState(const std::vector<std::string>& fluid_names) = 0;
1594 };
1595 
1599 void register_backend(const backend_families& bf, shared_ptr<AbstractStateGenerator> gen);
1600 
1601 template <class T>
1603 {
1604  public:
1606  register_backend(bf, shared_ptr<AbstractStateGenerator>(new T()));
1607  };
1608 };
1609 
1610 } /* namespace CoolProp */
1611 #endif /* ABSTRACTSTATE_H_ */