CoolProp  6.6.1dev
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_neff(void) {
164  throw NotImplementedError("calc_neff is not implemented for this backend");
165  };
167  virtual CoolPropDbl calc_umolar(void) {
168  throw NotImplementedError("calc_umolar is not implemented for this backend");
169  };
171  virtual CoolPropDbl calc_cpmolar(void) {
172  throw NotImplementedError("calc_cpmolar is not implemented for this backend");
173  };
176  throw NotImplementedError("calc_cpmolar_idealgas is not implemented for this backend");
177  };
179  virtual CoolPropDbl calc_cvmolar(void) {
180  throw NotImplementedError("calc_cvmolar is not implemented for this backend");
181  };
183  virtual CoolPropDbl calc_gibbsmolar(void) {
184  throw NotImplementedError("calc_gibbsmolar is not implemented for this backend");
185  };
188  throw NotImplementedError("calc_gibbsmolar_residual is not implemented for this backend");
189  };
192  throw NotImplementedError("calc_helmholtzmolar is not implemented for this backend");
193  };
196  throw NotImplementedError("calc_speed_sound is not implemented for this backend");
197  };
200  throw NotImplementedError("calc_isothermal_compressibility is not implemented for this backend");
201  };
204  throw NotImplementedError("calc_isobaric_expansion_coefficient is not implemented for this backend");
205  };
208  throw NotImplementedError("calc_isentropic_expansion_coefficient is not implemented for this backend");
209  };
211  virtual CoolPropDbl calc_viscosity(void) {
212  throw NotImplementedError("calc_viscosity is not implemented for this backend");
213  };
216  throw NotImplementedError("calc_conductivity is not implemented for this backend");
217  };
220  throw NotImplementedError("calc_surface_tension is not implemented for this backend");
221  };
223  virtual CoolPropDbl calc_molar_mass(void) {
224  throw NotImplementedError("calc_molar_mass is not implemented for this backend");
225  };
228  throw NotImplementedError("calc_acentric_factor is not implemented for this backend");
229  };
231  virtual CoolPropDbl calc_pressure(void) {
232  throw NotImplementedError("calc_pressure is not implemented for this backend");
233  };
236  throw NotImplementedError("calc_gas_constant is not implemented for this backend");
237  };
239  virtual CoolPropDbl calc_fugacity_coefficient(std::size_t i) {
240  throw NotImplementedError("calc_fugacity_coefficient is not implemented for this backend");
241  };
243  virtual std::vector<CoolPropDbl> calc_fugacity_coefficients() {
244  throw NotImplementedError("calc_fugacity_coefficients is not implemented for this backend");
245  };
247  virtual CoolPropDbl calc_fugacity(std::size_t i) {
248  throw NotImplementedError("calc_fugacity is not implemented for this backend");
249  };
251  virtual CoolPropDbl calc_chemical_potential(std::size_t i) {
252  throw NotImplementedError("calc_chemical_potential is not implemented for this backend");
253  };
255  virtual CoolPropDbl calc_PIP(void) {
256  throw NotImplementedError("calc_PIP is not implemented for this backend");
257  };
258 
259  // Excess properties
261  virtual void calc_excess_properties(void) {
262  throw NotImplementedError("calc_excess_properties is not implemented for this backend");
263  };
264 
265  // Derivatives of residual helmholtz energy
267  virtual CoolPropDbl calc_alphar(void) {
268  throw NotImplementedError("calc_alphar is not implemented for this backend");
269  };
272  throw NotImplementedError("calc_dalphar_dDelta is not implemented for this backend");
273  };
276  throw NotImplementedError("calc_dalphar_dTau is not implemented for this backend");
277  };
280  throw NotImplementedError("calc_d2alphar_dDelta2 is not implemented for this backend");
281  };
284  throw NotImplementedError("calc_d2alphar_dDelta_dTau is not implemented for this backend");
285  };
288  throw NotImplementedError("calc_d2alphar_dTau2 is not implemented for this backend");
289  };
292  throw NotImplementedError("calc_d3alphar_dDelta3 is not implemented for this backend");
293  };
296  throw NotImplementedError("calc_d3alphar_dDelta2_dTau is not implemented for this backend");
297  };
300  throw NotImplementedError("calc_d3alphar_dDelta_dTau2 is not implemented for this backend");
301  };
304  throw NotImplementedError("calc_d3alphar_dTau3 is not implemented for this backend");
305  };
306 
309  throw NotImplementedError("calc_d4alphar_dDelta4 is not implemented for this backend");
310  };
313  throw NotImplementedError("calc_d4alphar_dDelta3_dTau is not implemented for this backend");
314  };
317  throw NotImplementedError("calc_d4alphar_dDelta2_dTau2 is not implemented for this backend");
318  };
321  throw NotImplementedError("calc_d4alphar_dDelta_dTau3 is not implemented for this backend");
322  };
325  throw NotImplementedError("calc_d4alphar_dTau4 is not implemented for this backend");
326  };
327 
328  // Derivatives of ideal-gas helmholtz energy
330  virtual CoolPropDbl calc_alpha0(void) {
331  throw NotImplementedError("calc_alpha0 is not implemented for this backend");
332  };
335  throw NotImplementedError("calc_dalpha0_dDelta is not implemented for this backend");
336  };
339  throw NotImplementedError("calc_dalpha0_dTau is not implemented for this backend");
340  };
343  throw NotImplementedError("calc_d2alpha0_dDelta_dTau is not implemented for this backend");
344  };
347  throw NotImplementedError("calc_d2alpha0_dDelta2 is not implemented for this backend");
348  };
351  throw NotImplementedError("calc_d2alpha0_dTau2 is not implemented for this backend");
352  };
355  throw NotImplementedError("calc_d3alpha0_dDelta3 is not implemented for this backend");
356  };
359  throw NotImplementedError("calc_d3alpha0_dDelta2_dTau is not implemented for this backend");
360  };
363  throw NotImplementedError("calc_d3alpha0_dDelta_dTau2 is not implemented for this backend");
364  };
367  throw NotImplementedError("calc_d3alpha0_dTau3 is not implemented for this backend");
368  };
369 
370  virtual void calc_reducing_state(void) {
371  throw NotImplementedError("calc_reducing_state is not implemented for this backend");
372  };
373 
375  virtual CoolPropDbl calc_Tmax(void) {
376  throw NotImplementedError("calc_Tmax is not implemented for this backend");
377  };
379  virtual CoolPropDbl calc_Tmin(void) {
380  throw NotImplementedError("calc_Tmin is not implemented for this backend");
381  };
383  virtual CoolPropDbl calc_pmax(void) {
384  throw NotImplementedError("calc_pmax is not implemented for this backend");
385  };
386 
388  virtual CoolPropDbl calc_GWP20(void) {
389  throw NotImplementedError("calc_GWP20 is not implemented for this backend");
390  };
392  virtual CoolPropDbl calc_GWP100(void) {
393  throw NotImplementedError("calc_GWP100 is not implemented for this backend");
394  };
396  virtual CoolPropDbl calc_GWP500(void) {
397  throw NotImplementedError("calc_GWP500 is not implemented for this backend");
398  };
400  virtual CoolPropDbl calc_ODP(void) {
401  throw NotImplementedError("calc_ODP is not implemented for this backend");
402  };
405  throw NotImplementedError("calc_flame_hazard is not implemented for this backend");
406  };
409  throw NotImplementedError("calc_health_hazard is not implemented for this backend");
410  };
413  throw NotImplementedError("calc_physical_hazard is not implemented for this backend");
414  };
417  throw NotImplementedError("calc_dipole_moment is not implemented for this backend");
418  };
419 
423  virtual CoolPropDbl calc_second_partial_deriv(parameters Of1, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2);
424 
427  throw NotImplementedError("calc_reduced_density is not implemented for this backend");
428  };
431  throw NotImplementedError("calc_reciprocal_reduced_temperature is not implemented for this backend");
432  };
433 
435  virtual CoolPropDbl calc_Bvirial(void) {
436  throw NotImplementedError("calc_Bvirial is not implemented for this backend");
437  };
439  virtual CoolPropDbl calc_Cvirial(void) {
440  throw NotImplementedError("calc_Cvirial is not implemented for this backend");
441  };
444  throw NotImplementedError("calc_dBvirial_dT is not implemented for this backend");
445  };
448  throw NotImplementedError("calc_dCvirial_dT is not implemented for this backend");
449  };
452  throw NotImplementedError("calc_compressibility_factor is not implemented for this backend");
453  };
454 
456  virtual std::string calc_name(void) {
457  throw NotImplementedError("calc_name is not implemented for this backend");
458  };
460  virtual std::string calc_description(void) {
461  throw NotImplementedError("calc_description is not implemented for this backend");
462  };
463 
465  virtual CoolPropDbl calc_Ttriple(void) {
466  throw NotImplementedError("calc_Ttriple is not implemented for this backend");
467  };
469  virtual CoolPropDbl calc_p_triple(void) {
470  throw NotImplementedError("calc_p_triple is not implemented for this backend");
471  };
472 
474  virtual CoolPropDbl calc_T_critical(void) {
475  throw NotImplementedError("calc_T_critical is not implemented for this backend");
476  };
478  virtual CoolPropDbl calc_T_reducing(void) {
479  throw NotImplementedError("calc_T_reducing is not implemented for this backend");
480  };
482  virtual CoolPropDbl calc_p_critical(void) {
483  throw NotImplementedError("calc_p_critical is not implemented for this backend");
484  };
486  virtual CoolPropDbl calc_p_reducing(void) {
487  throw NotImplementedError("calc_p_reducing is not implemented for this backend");
488  };
491  throw NotImplementedError("calc_rhomolar_critical is not implemented for this backend");
492  };
495  throw NotImplementedError("calc_rhomass_critical is not implemented for this backend");
496  };
499  throw NotImplementedError("calc_rhomolar_reducing is not implemented for this backend");
500  };
502  virtual void calc_phase_envelope(const std::string& type) {
503  throw NotImplementedError("calc_phase_envelope is not implemented for this backend");
504  };
506  virtual CoolPropDbl calc_rhomass(void) {
507  return rhomolar() * molar_mass();
508  }
509  virtual CoolPropDbl calc_hmass(void) {
510  return hmolar() / molar_mass();
511  }
513  return hmolar_excess() / molar_mass();
514  }
515  virtual CoolPropDbl calc_smass(void) {
516  return smolar() / molar_mass();
517  }
519  return smolar_excess() / molar_mass();
520  }
521  virtual CoolPropDbl calc_cpmass(void) {
522  return cpmolar() / molar_mass();
523  }
524  virtual CoolPropDbl calc_cp0mass(void) {
525  return cp0molar() / molar_mass();
526  }
527  virtual CoolPropDbl calc_cvmass(void) {
528  return cvmolar() / molar_mass();
529  }
530  virtual CoolPropDbl calc_umass(void) {
531  return umolar() / molar_mass();
532  }
534  return umolar_excess() / molar_mass();
535  }
536  virtual CoolPropDbl calc_gibbsmass(void) {
537  return gibbsmolar() / molar_mass();
538  }
540  return gibbsmolar_excess() / molar_mass();
541  }
543  return helmholtzmolar() / molar_mass();
544  }
546  return helmholtzmolar_excess() / molar_mass();
547  }
549  return volumemolar_excess() / molar_mass();
550  }
551 
553  virtual void update_states(void) {
554  throw NotImplementedError("This backend does not implement update_states function");
555  };
556 
557  virtual CoolPropDbl calc_melting_line(int param, int given, CoolPropDbl value) {
558  throw NotImplementedError("This backend does not implement calc_melting_line function");
559  };
560 
565  virtual CoolPropDbl calc_saturation_ancillary(parameters param, int Q, parameters given, double value) {
566  throw NotImplementedError("This backend does not implement calc_saturation_ancillary");
567  };
568 
570  virtual phases calc_phase(void) {
571  throw NotImplementedError("This backend does not implement calc_phase function");
572  };
575  throw NotImplementedError("This backend does not implement calc_specify_phase function");
576  };
578  virtual void calc_unspecify_phase(void) {
579  throw NotImplementedError("This backend does not implement calc_unspecify_phase function");
580  };
582  virtual std::vector<std::string> calc_fluid_names(void) {
583  throw NotImplementedError("This backend does not implement calc_fluid_names function");
584  };
587  virtual const CoolProp::SimpleState& calc_state(const std::string& state) {
588  throw NotImplementedError("calc_state is not implemented for this backend");
589  };
590 
592  throw NotImplementedError("calc_phase_envelope_data is not implemented for this backend");
593  };
594 
595  virtual std::vector<CoolPropDbl> calc_mole_fractions_liquid(void) {
596  throw NotImplementedError("calc_mole_fractions_liquid is not implemented for this backend");
597  };
598  virtual std::vector<CoolPropDbl> calc_mole_fractions_vapor(void) {
599  throw NotImplementedError("calc_mole_fractions_vapor is not implemented for this backend");
600  };
601  virtual const std::vector<CoolPropDbl> calc_mass_fractions(void) {
602  throw NotImplementedError("calc_mass_fractions is not implemented for this backend");
603  };
604 
607  throw NotImplementedError("calc_fraction_min is not implemented for this backend");
608  };
611  throw NotImplementedError("calc_fraction_max is not implemented for this backend");
612  };
613  virtual CoolPropDbl calc_T_freeze(void) {
614  throw NotImplementedError("calc_T_freeze is not implemented for this backend");
615  };
616 
618  throw NotImplementedError("calc_first_saturation_deriv is not implemented for this backend");
619  };
621  throw NotImplementedError("calc_second_saturation_deriv is not implemented for this backend");
622  };
624  throw NotImplementedError("calc_first_two_phase_deriv is not implemented for this backend");
625  };
627  throw NotImplementedError("calc_second_two_phase_deriv is not implemented for this backend");
628  };
630  throw NotImplementedError("calc_first_two_phase_deriv_splined is not implemented for this backend");
631  };
632 
634  throw NotImplementedError("calc_saturated_liquid_keyed_output is not implemented for this backend");
635  };
637  throw NotImplementedError("calc_saturated_vapor_keyed_output is not implemented for this backend");
638  };
639  virtual void calc_ideal_curve(const std::string& type, std::vector<double>& T, std::vector<double>& p) {
640  throw NotImplementedError("calc_ideal_curve is not implemented for this backend");
641  };
642 
644  virtual CoolPropDbl calc_T(void) {
645  return _T;
646  }
648  virtual CoolPropDbl calc_rhomolar(void) {
649  return _rhomolar;
650  }
651 
653  virtual double calc_tangent_plane_distance(const double T, const double p, const std::vector<double>& w, const double rhomolar_guess) {
654  throw NotImplementedError("calc_tangent_plane_distance is not implemented for this backend");
655  };
656 
658  virtual void calc_true_critical_point(double& T, double& rho) {
659  throw NotImplementedError("calc_true_critical_point is not implemented for this backend");
660  };
661 
662  virtual void calc_conformal_state(const std::string& reference_fluid, CoolPropDbl& T, CoolPropDbl& rhomolar) {
663  throw NotImplementedError("calc_conformal_state is not implemented for this backend");
664  };
665 
666  virtual void calc_viscosity_contributions(CoolPropDbl& dilute, CoolPropDbl& initial_density, CoolPropDbl& residual, CoolPropDbl& critical) {
667  throw NotImplementedError("calc_viscosity_contributions is not implemented for this backend");
668  };
669  virtual void calc_conductivity_contributions(CoolPropDbl& dilute, CoolPropDbl& initial_density, CoolPropDbl& residual, CoolPropDbl& critical) {
670  throw NotImplementedError("calc_conductivity_contributions is not implemented for this backend");
671  };
672  virtual std::vector<CriticalState> calc_all_critical_points(void) {
673  throw NotImplementedError("calc_all_critical_points is not implemented for this backend");
674  };
675  virtual void calc_build_spinodal() {
676  throw NotImplementedError("calc_build_spinodal is not implemented for this backend");
677  };
679  throw NotImplementedError("calc_get_spinodal_data is not implemented for this backend");
680  };
681  virtual void calc_criticality_contour_values(double& L1star, double& M1star) {
682  throw NotImplementedError("calc_criticality_contour_values is not implemented for this backend");
683  };
684 
686  virtual void mass_to_molar_inputs(CoolProp::input_pairs& input_pair, CoolPropDbl& value1, CoolPropDbl& value2);
687 
689  virtual void calc_change_EOS(const std::size_t i, const std::string& EOS_name) {
690  throw NotImplementedError("calc_change_EOS is not implemented for this backend");
691  };
692 
693  public:
695  clear();
696  }
697  virtual ~AbstractState(){};
698 
700 
706  static AbstractState* factory(const std::string& backend, const std::string& fluid_names) {
707  return factory(backend, strsplit(fluid_names, '&'));
708  };
709 
733  static AbstractState* factory(const std::string& backend, const std::vector<std::string>& fluid_names);
734 
737  _T = T;
738  }
739 
744  virtual std::string backend_name(void) = 0;
745 
746  // The derived classes must implement this function to define whether they use mole fractions (true) or mass fractions (false)
747  virtual bool using_mole_fractions(void) = 0;
748  virtual bool using_mass_fractions(void) = 0;
749  virtual bool using_volu_fractions(void) = 0;
750 
751  virtual void set_mole_fractions(const std::vector<CoolPropDbl>& mole_fractions) = 0;
752  virtual void set_mass_fractions(const std::vector<CoolPropDbl>& mass_fractions) = 0;
753  virtual void set_volu_fractions(const std::vector<CoolPropDbl>& mass_fractions) {
754  throw NotImplementedError("Volume composition has not been implemented.");
755  }
756 
776  virtual void set_reference_stateS(const std::string& reference_state) {
777  throw NotImplementedError(
778  "Setting reference state has not been implemented for this backend. Try using CoolProp::set_reference_stateD instead.");
779  }
780 
786  virtual void set_reference_stateD(double T, double rhomolar, double hmolar0, double smolar0) {
787  throw NotImplementedError(
788  "Setting reference state has not been implemented for this backend. Try using CoolProp::set_reference_stateD instead.");
789  }
790 
791 #ifndef COOLPROPDBL_MAPS_TO_DOUBLE
792  void set_mole_fractions(const std::vector<double>& mole_fractions) {
793  set_mole_fractions(std::vector<CoolPropDbl>(mole_fractions.begin(), mole_fractions.end()));
794  };
795  void set_mass_fractions(const std::vector<double>& mass_fractions) {
796  set_mass_fractions(std::vector<CoolPropDbl>(mass_fractions.begin(), mass_fractions.end()));
797  };
798  void set_volu_fractions(const std::vector<double>& volu_fractions) {
799  set_volu_fractions(std::vector<CoolPropDbl>(volu_fractions.begin(), volu_fractions.end()));
800  };
801 #endif
802 
803 #ifdef EMSCRIPTEN
804  void set_mole_fractions_double(const std::vector<double>& mole_fractions) {
805  set_mole_fractions(std::vector<CoolPropDbl>(mole_fractions.begin(), mole_fractions.end()));
806  };
807 #endif
808 
810  std::vector<CoolPropDbl> mole_fractions_liquid(void) {
812  };
814  std::vector<double> mole_fractions_liquid_double(void) {
815  std::vector<CoolPropDbl> x = calc_mole_fractions_liquid();
816  return std::vector<double>(x.begin(), x.end());
817  };
818 
820  std::vector<CoolPropDbl> mole_fractions_vapor(void) {
821  return calc_mole_fractions_vapor();
822  };
824  std::vector<double> mole_fractions_vapor_double(void) {
825  std::vector<CoolPropDbl> y = calc_mole_fractions_vapor();
826  return std::vector<double>(y.begin(), y.end());
827  };
828 
830  virtual const std::vector<CoolPropDbl>& get_mole_fractions(void) = 0;
832  virtual const std::vector<CoolPropDbl> get_mass_fractions(void) {
833  return this->calc_mass_fractions();
834  };
835 
837  virtual void update(CoolProp::input_pairs input_pair, double Value1, double Value2) = 0;
838 
841  virtual void update_with_guesses(CoolProp::input_pairs input_pair, double Value1, double Value2, const GuessesStructure& guesses) {
842  throw NotImplementedError("update_with_guesses is not implemented for this backend");
843  };
844 
848  virtual bool available_in_high_level(void) {
849  return true;
850  }
851 
853  virtual std::string fluid_param_string(const std::string&) {
854  throw NotImplementedError("fluid_param_string has not been implemented for this backend");
855  }
856 
858  std::vector<std::string> fluid_names(void);
859 
864  virtual const double get_fluid_constant(std::size_t i, parameters param) const {
865  throw NotImplementedError("get_fluid_constant is not implemented for this backend");
866  };
867  ;
868 
870  virtual void set_binary_interaction_double(const std::string& CAS1, const std::string& CAS2, 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_double(const std::size_t i, const std::size_t j, const std::string& parameter, const double value) {
875  throw NotImplementedError("set_binary_interaction_double is not implemented for this backend");
876  };
878  virtual void set_binary_interaction_string(const std::string& CAS1, const std::string& CAS2, const std::string& parameter,
879  const std::string& value) {
880  throw NotImplementedError("set_binary_interaction_string is not implemented for this backend");
881  };
883  virtual void set_binary_interaction_string(const std::size_t i, const std::size_t j, const std::string& parameter, const std::string& value) {
884  throw NotImplementedError("set_binary_interaction_string is not implemented for this backend");
885  };
887  virtual double get_binary_interaction_double(const std::string& CAS1, const std::string& CAS2, const std::string& parameter) {
888  throw NotImplementedError("get_binary_interaction_double is not implemented for this backend");
889  };
891  virtual double get_binary_interaction_double(const std::size_t i, const std::size_t j, const std::string& parameter) {
892  throw NotImplementedError("get_binary_interaction_double is not implemented for this backend");
893  };
895  virtual std::string get_binary_interaction_string(const std::string& CAS1, const std::string& CAS2, const std::string& parameter) {
896  throw NotImplementedError("get_binary_interaction_string is not implemented for this backend");
897  };
899  virtual void apply_simple_mixing_rule(std::size_t i, std::size_t j, const std::string& model) {
900  throw NotImplementedError("apply_simple_mixing_rule is not implemented for this backend");
901  };
903  virtual void set_cubic_alpha_C(const size_t i, const std::string& parameter, const double c1, const double c2, const double c3) {
904  throw ValueError("set_cubic_alpha_C only defined for cubic backends");
905  };
907  virtual void set_fluid_parameter_double(const size_t i, const std::string& parameter, const double value) {
908  throw ValueError("set_fluid_parameter_double only defined for cubic backends");
909  };
911  virtual double get_fluid_parameter_double(const size_t i, const std::string& parameter) {
912  throw ValueError("get_fluid_parameter_double only defined for cubic backends");
913  };
914 
916  virtual bool clear();
918  virtual bool clear_comp_change();
919 
924  return _reducing;
925  };
926 
928  const CoolProp::SimpleState& get_state(const std::string& state) {
929  return calc_state(state);
930  };
931 
933  double Tmin(void);
935  double Tmax(void);
937  double pmax(void);
939  double Ttriple(void);
940 
942  phases phase(void) {
943  return calc_phase();
944  };
948  };
950  void unspecify_phase(void) {
952  };
953 
955  double T_critical(void);
957  double p_critical(void);
959  double rhomolar_critical(void);
961  double rhomass_critical(void);
962 
964  std::vector<CriticalState> all_critical_points(void) {
965  return calc_all_critical_points();
966  };
967 
969  void build_spinodal() {
971  };
972 
975  return calc_get_spinodal_data();
976  };
977 
979  void criticality_contour_values(double& L1star, double& M1star) {
980  return calc_criticality_contour_values(L1star, M1star);
981  }
982 
1006  double tangent_plane_distance(const double T, const double p, const std::vector<double>& w, const double rhomolar_guess = -1) {
1007  return calc_tangent_plane_distance(T, p, w, rhomolar_guess);
1008  };
1009 
1011  double T_reducing(void);
1013  double rhomolar_reducing(void);
1015  double rhomass_reducing(void);
1016 
1018  double p_triple(void);
1019 
1021  std::string name() {
1022  return calc_name();
1023  };
1025  std::string description() {
1026  return calc_description();
1027  };
1028 
1030  double dipole_moment() {
1031  return calc_dipole_moment();
1032  }
1033 
1034  // ----------------------------------------
1035  // Bulk properties - temperature and density are directly calculated every time
1036  // All other parameters are calculated on an as-needed basis
1037  // ----------------------------------------
1039  double keyed_output(parameters key);
1041  double trivial_keyed_output(parameters key);
1045  };
1049  };
1050 
1052  double T(void) {
1053  return calc_T();
1054  };
1056  double rhomolar(void) {
1057  return calc_rhomolar();
1058  };
1060  double rhomass(void) {
1061  return calc_rhomass();
1062  };
1064  double p(void) {
1065  return _p;
1066  };
1068  double Q(void) {
1069  return _Q;
1070  };
1072  double tau(void);
1074  double delta(void);
1076  double molar_mass(void);
1078  double acentric_factor(void);
1080  double gas_constant(void);
1082  double Bvirial(void);
1084  double dBvirial_dT(void);
1086  double Cvirial(void);
1088  double dCvirial_dT(void);
1090  double compressibility_factor(void);
1092  double hmolar(void);
1094  double hmolar_residual(void);
1096  double hmass(void) {
1097  return calc_hmass();
1098  };
1100  double hmolar_excess(void);
1102  double hmass_excess(void) {
1103  return calc_hmass_excess();
1104  };
1106  double smolar(void);
1108  double smolar_residual(void);
1110  double neff(void);
1112  double smass(void) {
1113  return calc_smass();
1114  };
1116  double smolar_excess(void);
1118  double smass_excess(void) {
1119  return calc_smass_excess();
1120  };
1122  double umolar(void);
1124  double umass(void) {
1125  return calc_umass();
1126  };
1128  double umolar_excess(void);
1130  double umass_excess(void) {
1131  return calc_umass_excess();
1132  };
1134  double cpmolar(void);
1136  double cpmass(void) {
1137  return calc_cpmass();
1138  };
1140  double cp0molar(void);
1142  double cp0mass(void) {
1143  return calc_cp0mass();
1144  };
1146  double cvmolar(void);
1148  double cvmass(void) {
1149  return calc_cvmass();
1150  };
1152  double gibbsmolar(void);
1154  double gibbsmolar_residual(void);
1156  double gibbsmass(void) {
1157  return calc_gibbsmass();
1158  };
1160  double gibbsmolar_excess(void);
1162  double gibbsmass_excess(void) {
1163  return calc_gibbsmass_excess();
1164  };
1166  double helmholtzmolar(void);
1168  double helmholtzmass(void) {
1169  return calc_helmholtzmass();
1170  };
1172  double helmholtzmolar_excess(void);
1174  double helmholtzmass_excess(void) {
1175  return calc_helmholtzmass_excess();
1176  };
1178  double volumemolar_excess(void);
1180  double volumemass_excess(void) {
1181  return calc_volumemass_excess();
1182  };
1184  double speed_sound(void);
1186  double isothermal_compressibility(void);
1188  double isobaric_expansion_coefficient(void);
1190  double isentropic_expansion_coefficient(void);
1192  double fugacity_coefficient(std::size_t i);
1194  std::vector<double> fugacity_coefficients();
1196  double fugacity(std::size_t i);
1198  double chemical_potential(std::size_t i);
1209  double PIP() {
1210  return calc_PIP();
1211  };
1212 
1214  void true_critical_point(double& T, double& rho) {
1216  }
1217 
1225  void ideal_curve(const std::string& type, std::vector<double>& T, std::vector<double>& p) {
1226  calc_ideal_curve(type, T, p);
1227  };
1228 
1229  // ----------------------------------------
1230  // Partial derivatives
1231  // ----------------------------------------
1232 
1238  return calc_first_partial_deriv(Of, Wrt, Constant);
1239  };
1240 
1265  return calc_second_partial_deriv(Of1, Wrt1, Constant1, Wrt2, Constant2);
1266  };
1267 
1290  return calc_first_saturation_deriv(Of1, Wrt1);
1291  };
1292 
1311  return calc_second_saturation_deriv(Of1, Wrt1, Wrt2);
1312  };
1313 
1334  return calc_first_two_phase_deriv(Of, Wrt, Constant);
1335  };
1336 
1353  double second_two_phase_deriv(parameters Of, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2) {
1354  return calc_second_two_phase_deriv(Of, Wrt1, Constant1, Wrt2, Constant2);
1355  };
1356 
1377  double first_two_phase_deriv_splined(parameters Of, parameters Wrt, parameters Constant, double x_end) {
1378  return calc_first_two_phase_deriv_splined(Of, Wrt, Constant, x_end);
1379  };
1380 
1381  // ----------------------------------------
1382  // Phase envelope for mixtures
1383  // ----------------------------------------
1384 
1390  void build_phase_envelope(const std::string& type = "");
1395  return calc_phase_envelope_data();
1396  };
1397 
1398  // ----------------------------------------
1399  // Ancillary equations
1400  // ----------------------------------------
1401 
1403  virtual bool has_melting_line(void) {
1404  return false;
1405  };
1410  double melting_line(int param, int given, double value);
1416  double saturation_ancillary(parameters param, int Q, parameters given, double value);
1417 
1418  // ----------------------------------------
1419  // Transport properties
1420  // ----------------------------------------
1422  double viscosity(void);
1424  void viscosity_contributions(CoolPropDbl& dilute, CoolPropDbl& initial_density, CoolPropDbl& residual, CoolPropDbl& critical) {
1425  calc_viscosity_contributions(dilute, initial_density, residual, critical);
1426  };
1428  double conductivity(void);
1430  void conductivity_contributions(CoolPropDbl& dilute, CoolPropDbl& initial_density, CoolPropDbl& residual, CoolPropDbl& critical) {
1431  calc_conductivity_contributions(dilute, initial_density, residual, critical);
1432  };
1434  double surface_tension(void);
1436  double Prandtl(void) {
1437  return cpmass() * viscosity() / conductivity();
1438  };
1445  void conformal_state(const std::string& reference_fluid, CoolPropDbl& T, CoolPropDbl& rhomolar) {
1446  return calc_conformal_state(reference_fluid, T, rhomolar);
1447  };
1448 
1453  void change_EOS(const std::size_t i, const std::string& EOS_name) {
1454  calc_change_EOS(i, EOS_name);
1455  }
1456 
1457  // ----------------------------------------
1458  // Helmholtz energy and derivatives
1459  // ----------------------------------------
1462  if (!_alpha0) _alpha0 = calc_alpha0();
1463  return _alpha0;
1464  };
1468  return _dalpha0_dDelta;
1469  };
1473  return _dalpha0_dTau;
1474  };
1478  return _d2alpha0_dDelta2;
1479  };
1483  return _d2alpha0_dDelta_dTau;
1484  };
1488  return _d2alpha0_dTau2;
1489  };
1493  return _d3alpha0_dTau3;
1494  };
1498  return _d3alpha0_dDelta_dTau2;
1499  };
1503  return _d3alpha0_dDelta2_dTau;
1504  };
1508  return _d3alpha0_dDelta3;
1509  };
1510 
1513  if (!_alphar) _alphar = calc_alphar();
1514  return _alphar;
1515  };
1519  return _dalphar_dDelta;
1520  };
1524  return _dalphar_dTau;
1525  };
1529  return _d2alphar_dDelta2;
1530  };
1534  return _d2alphar_dDelta_dTau;
1535  };
1539  return _d2alphar_dTau2;
1540  };
1544  return _d3alphar_dDelta3;
1545  };
1549  return _d3alphar_dDelta2_dTau;
1550  };
1554  return _d3alphar_dDelta_dTau2;
1555  };
1559  return _d3alphar_dTau3;
1560  };
1564  return _d4alphar_dDelta4;
1565  };
1569  return _d4alphar_dDelta3_dTau;
1570  };
1574  return _d4alphar_dDelta2_dTau2;
1575  };
1579  return _d4alphar_dDelta_dTau3;
1580  };
1584  return _d4alphar_dTau4;
1585  };
1586 };
1587 
1596 {
1597  public:
1598  virtual AbstractState* get_AbstractState(const std::vector<std::string>& fluid_names) = 0;
1600 };
1601 
1605 void register_backend(const backend_families& bf, shared_ptr<AbstractStateGenerator> gen);
1606 
1607 template <class T>
1609 {
1610  public:
1612  register_backend(bf, shared_ptr<AbstractStateGenerator>(new T()));
1613  };
1614 };
1615 
1616 } /* namespace CoolProp */
1617 #endif /* ABSTRACTSTATE_H_ */