CoolProp
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
CoolProp::AbstractState Class Referenceabstract

Detailed Description

The mother of all state classes.

This class provides the basic properties based on interrelations of the properties, their derivatives and the Helmholtz energy terms. It does not provide the mechanism to update the values. This has to be implemented in a subclass. Most functions are defined as virtual functions allowing us redefine them later, for example to implement the TTSE technique. The functions defined here are always used as a fall-back.

This base class does not perform any checks on the two-phase conditions and alike. Most of the functions defined here only apply to compressible single state substances. Make sure you are aware of all the assumptions we made when using this class.

Add build table function to Abstract State Interpolator inherit AS implemented by TTSE BICUBIC

#include <AbstractState.h>

Inheritance diagram for CoolProp::AbstractState:
CoolProp::HelmholtzEOSMixtureBackend CoolProp::IF97Backend CoolProp::IncompressibleBackend CoolProp::REFPROPMixtureBackend CoolProp::TabularBackend CoolProp::AbstractCubicBackend CoolProp::HelmholtzEOSBackend CoolProp::REFPROPBackend CoolProp::BicubicBackend CoolProp::TTSEBackend CoolProp::PengRobinsonBackend CoolProp::SRKBackend CoolProp::VTPRBackend

Public Member Functions

void set_T (CoolPropDbl T)
 Set the internal variable T without a flash call (expert use only!)
 
virtual std::string backend_name (void)=0
 Get a string representation of the backend - for instance "HelmholtzEOSMixtureBackend" for the core mixture model in CoolProp. More...
 
virtual bool using_mole_fractions (void)=0
 
virtual bool using_mass_fractions (void)=0
 
virtual bool using_volu_fractions (void)=0
 
virtual void set_mole_fractions (const std::vector< CoolPropDbl > &mole_fractions)=0
 
virtual void set_mass_fractions (const std::vector< CoolPropDbl > &mass_fractions)=0
 
virtual void set_volu_fractions (const std::vector< CoolPropDbl > &mass_fractions)
 
virtual void set_reference_stateS (const std::string &reference_state)
 Set the reference state based on a string representation. More...
 
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 density. More...
 
std::vector< CoolPropDbl > mole_fractions_liquid (void)
 Get the mole fractions of the equilibrium liquid phase.
 
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)
 
std::vector< CoolPropDbl > mole_fractions_vapor (void)
 Get the mole fractions of the equilibrium vapor phase.
 
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)
 
virtual const std::vector< CoolPropDbl > & get_mole_fractions (void)=0
 Get the mole fractions of the fluid.
 
virtual const std::vector< CoolPropDbl > get_mass_fractions (void)
 Get the mass fractions of the fluid.
 
virtual void update (CoolProp::input_pairs input_pair, double Value1, double Value2)=0
 Update the state using two state variables.
 
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 be used - this is backend dependent.
 
virtual bool available_in_high_level (void)
 A function that says whether the backend instance can be instantiated in the high-level interface In general this should be true, except for some other backends (especially the tabular backends) To disable use in high-level interface, implement this function and return false.
 
virtual std::string fluid_param_string (const std::string &)
 Return a string from the backend for the mixture/fluid - backend dependent - could be CAS #, name, etc.
 
std::vector< std::string > fluid_names (void)
 Return a vector of strings of the fluid names that are in use.
 
virtual const double get_fluid_constant (std::size_t i, parameters param) const
 Get a constant for one of the fluids forming this mixture. More...
 
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!!!)
 
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!!!)
 
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!!!)
 
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!!!)
 
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!!!)
 
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!!!)
 
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!!!)
 
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!!!)
 
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's constants:
 
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)
 
virtual double get_fluid_parameter_double (const size_t i, const std::string &parameter)
 Double fluid parameter (currently the volume translation parameter for cubic)
 
virtual bool clear ()
 Clear all the cached values. More...
 
virtual const CoolProp::SimpleStateget_reducing_state ()
 Get the state that is used in the equation of state or mixture model to reduce the state. More...
 
const CoolProp::SimpleStateget_state (const std::string &state)
 Get a desired state point - backend dependent.
 
double Tmin (void)
 Get the minimum temperature in K.
 
double Tmax (void)
 Get the maximum temperature in K.
 
double pmax (void)
 Get the maximum pressure in Pa.
 
double Ttriple (void)
 Get the triple point temperature in K.
 
phases phase (void)
 Get the phase of the state.
 
void specify_phase (phases phase)
 Specify the phase for all further calculations with this state class.
 
void unspecify_phase (void)
 Unspecify the phase and go back to calculating it based on the inputs.
 
double T_critical (void)
 Return the critical temperature in K.
 
double p_critical (void)
 Return the critical pressure in Pa.
 
double rhomolar_critical (void)
 Return the critical molar density in mol/m^3.
 
double rhomass_critical (void)
 Return the critical mass density in kg/m^3.
 
std::vector< CriticalStateall_critical_points (void)
 Return the vector of critical points, including points that are unstable or correspond to negative pressure.
 
void build_spinodal ()
 Construct the spinodal curve for the mixture (or pure fluid)
 
SpinodalData get_spinodal_data ()
 Get the data from the spinodal curve constructed in the call to build_spinodal()
 
void criticality_contour_values (double &L1star, double &M1star)
 Calculate the criticality contour values \(\mathcal{L}_1^*\) and \(\mathcal{M}_1^*\).
 
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. More...
 
double T_reducing (void)
 Return the reducing point temperature in K.
 
double rhomolar_reducing (void)
 Return the molar density at the reducing point in mol/m^3.
 
double rhomass_reducing (void)
 Return the mass density at the reducing point in kg/m^3.
 
double p_triple (void)
 Return the triple point pressure in Pa.
 
std::string name ()
 Return the name - backend dependent.
 
std::string description ()
 Return the description - backend dependent.
 
double dipole_moment ()
 Return the dipole moment in C-m (1 D = 3.33564e-30 C-m)
 
double keyed_output (parameters key)
 Retrieve a value by key.
 
double trivial_keyed_output (parameters key)
 A trivial keyed output like molar mass that does not depend on the state.
 
double saturated_liquid_keyed_output (parameters key)
 Get an output from the saturated liquid state by key.
 
double saturated_vapor_keyed_output (parameters key)
 Get an output from the saturated vapor state by key.
 
double T (void)
 Return the temperature in K.
 
double rhomolar (void)
 Return the molar density in mol/m^3.
 
double rhomass (void)
 Return the mass density in kg/m^3.
 
double p (void)
 Return the pressure in Pa.
 
double Q (void)
 Return the vapor quality (mol/mol); Q = 0 for saturated liquid.
 
double tau (void)
 Return the reciprocal of the reduced temperature ( \(\tau = T_c/T\))
 
double delta (void)
 Return the reduced density ( \(\delta = \rho/\rho_c\))
 
double molar_mass (void)
 Return the molar mass in kg/mol.
 
double acentric_factor (void)
 Return the acentric factor.
 
double gas_constant (void)
 Return the mole-fraction weighted gas constant in J/mol/K.
 
double Bvirial (void)
 Return the B virial coefficient.
 
double dBvirial_dT (void)
 Return the derivative of the B virial coefficient with respect to temperature.
 
double Cvirial (void)
 Return the C virial coefficient.
 
double dCvirial_dT (void)
 Return the derivative of the C virial coefficient with respect to temperature.
 
double compressibility_factor (void)
 Return the compressibility factor \( Z = p/(rho R T) \).
 
double hmolar (void)
 Return the molar enthalpy in J/mol.
 
double hmass (void)
 Return the mass enthalpy in J/kg.
 
double hmolar_excess (void)
 Return the excess molar enthalpy in J/mol.
 
double hmass_excess (void)
 Return the excess mass enthalpy in J/kg.
 
double smolar (void)
 Return the molar entropy in J/mol/K.
 
double smass (void)
 Return the molar entropy in J/kg/K.
 
double smolar_excess (void)
 Return the molar entropy in J/mol/K.
 
double smass_excess (void)
 Return the molar entropy in J/kg/K.
 
double umolar (void)
 Return the molar internal energy in J/mol.
 
double umass (void)
 Return the mass internal energy in J/kg.
 
double umolar_excess (void)
 Return the excess internal energy in J/mol.
 
double umass_excess (void)
 Return the excess internal energy in J/kg.
 
double cpmolar (void)
 Return the molar constant pressure specific heat in J/mol/K.
 
double cpmass (void)
 Return the mass constant pressure specific heat in J/kg/K.
 
double cp0molar (void)
 Return the molar constant pressure specific heat for ideal gas part only in J/mol/K.
 
double cp0mass (void)
 Return the mass constant pressure specific heat for ideal gas part only in J/kg/K.
 
double cvmolar (void)
 Return the molar constant volume specific heat in J/mol/K.
 
double cvmass (void)
 Return the mass constant volume specific heat in J/kg/K.
 
double gibbsmolar (void)
 Return the Gibbs energy in J/mol.
 
double gibbsmass (void)
 Return the Gibbs energy in J/kg.
 
double gibbsmolar_excess (void)
 Return the excess Gibbs energy in J/mol.
 
double gibbsmass_excess (void)
 Return the excess Gibbs energy in J/kg.
 
double helmholtzmolar (void)
 Return the Helmholtz energy in J/mol.
 
double helmholtzmass (void)
 Return the Helmholtz energy in J/kg.
 
double helmholtzmolar_excess (void)
 Return the excess Helmholtz energy in J/mol.
 
double helmholtzmass_excess (void)
 Return the excess Helmholtz energy in J/kg.
 
double volumemolar_excess (void)
 Return the excess volume in m^3/mol.
 
double volumemass_excess (void)
 Return the excess volume in m^3/kg.
 
double speed_sound (void)
 Return the speed of sound in m/s.
 
double isothermal_compressibility (void)
 Return the isothermal compressibility \( \kappa = -\frac{1}{v}\left.\frac{\partial v}{\partial p}\right|_T=\frac{1}{\rho}\left.\frac{\partial \rho}{\partial p}\right|_T\) in 1/Pa.
 
double isobaric_expansion_coefficient (void)
 Return the isobaric expansion coefficient \( \beta = \frac{1}{v}\left.\frac{\partial v}{\partial T}\right|_p = -\frac{1}{\rho}\left.\frac{\partial \rho}{\partial T}\right|_p\) in 1/K.
 
double isentropic_expansion_coefficient (void)
 Return the isentropic expansion coefficient \( \kappa_s = -\frac{c_p}{c_v}\frac{v}{p}\left.\frac{\partial p}{\partial v}\right|_T = \frac{\rho}{p}\left.\frac{\partial p}{\partial \rho}\right|_s\).
 
double fugacity_coefficient (std::size_t i)
 Return the fugacity coefficient of the i-th component of the mixture.
 
double fugacity (std::size_t i)
 Return the fugacity of the i-th component of the mixture.
 
double chemical_potential (std::size_t i)
 Return the chemical potential of the i-th component of the mixture.
 
double fundamental_derivative_of_gas_dynamics (void)
 Return the fundamental derivative of gas dynamics \( \Gamma \). More...
 
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".
 
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.
 
void ideal_curve (const std::string &type, std::vector< double > &T, std::vector< double > &p)
 Calculate an ideal curve for a pure fluid. More...
 
CoolPropDbl first_partial_deriv (parameters Of, parameters Wrt, parameters Constant)
 The first partial derivative in homogeneous phases. More...
 
CoolPropDbl second_partial_deriv (parameters Of1, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2)
 The second partial derivative in homogeneous phases. More...
 
CoolPropDbl first_saturation_deriv (parameters Of1, parameters Wrt1)
 The first partial derivative along the saturation curve. More...
 
CoolPropDbl second_saturation_deriv (parameters Of1, parameters Wrt1, parameters Wrt2)
 The second partial derivative along the saturation curve. More...
 
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. More...
 
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. More...
 
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. More...
 
void build_phase_envelope (const std::string &type="")
 Construct the phase envelope for a mixture. More...
 
const CoolProp::PhaseEnvelopeDataget_phase_envelope_data ()
 After having calculated the phase envelope, return the phase envelope data.
 
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 class.
 
double melting_line (int param, int given, double value)
 Return a value from the melting line. More...
 
double saturation_ancillary (parameters param, int Q, parameters given, double value)
 Return the value from a saturation ancillary curve (if the backend implements it) More...
 
double viscosity (void)
 Return the viscosity in Pa-s.
 
void viscosity_contributions (CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical)
 Return the viscosity contributions, each in Pa-s.
 
double conductivity (void)
 Return the thermal conductivity in W/m/K.
 
void conductivity_contributions (CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical)
 Return the thermal conductivity contributions, each in W/m/K.
 
double surface_tension (void)
 Return the surface tension in N/m.
 
double Prandtl (void)
 Return the Prandtl number (dimensionless)
 
void conformal_state (const std::string &reference_fluid, CoolPropDbl &T, CoolPropDbl &rhomolar)
 Find the conformal state needed for ECS. More...
 
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. More...
 
CoolPropDbl alpha0 (void)
 Return the term \( \alpha^0 \).
 
CoolPropDbl dalpha0_dDelta (void)
 Return the term \( \alpha^0_{\delta} \).
 
CoolPropDbl dalpha0_dTau (void)
 Return the term \( \alpha^0_{\tau} \).
 
CoolPropDbl d2alpha0_dDelta2 (void)
 Return the term \( \alpha^0_{\delta\delta} \).
 
CoolPropDbl d2alpha0_dDelta_dTau (void)
 Return the term \( \alpha^0_{\delta\tau} \).
 
CoolPropDbl d2alpha0_dTau2 (void)
 Return the term \( \alpha^0_{\tau\tau} \).
 
CoolPropDbl d3alpha0_dTau3 (void)
 Return the term \( \alpha^0_{\tau\tau\tau} \).
 
CoolPropDbl d3alpha0_dDelta_dTau2 (void)
 Return the term \( \alpha^0_{\delta\tau\tau} \).
 
CoolPropDbl d3alpha0_dDelta2_dTau (void)
 Return the term \( \alpha^0_{\delta\delta\tau} \).
 
CoolPropDbl d3alpha0_dDelta3 (void)
 Return the term \( \alpha^0_{\delta\delta\delta} \).
 
CoolPropDbl alphar (void)
 Return the term \( \alpha^r \).
 
CoolPropDbl dalphar_dDelta (void)
 Return the term \( \alpha^r_{\delta} \).
 
CoolPropDbl dalphar_dTau (void)
 Return the term \( \alpha^r_{\tau} \).
 
CoolPropDbl d2alphar_dDelta2 (void)
 Return the term \( \alpha^r_{\delta\delta} \).
 
CoolPropDbl d2alphar_dDelta_dTau (void)
 Return the term \( \alpha^r_{\delta\tau} \).
 
CoolPropDbl d2alphar_dTau2 (void)
 Return the term \( \alpha^r_{\tau\tau} \).
 
CoolPropDbl d3alphar_dDelta3 (void)
 Return the term \( \alpha^r_{\delta\delta\delta} \).
 
CoolPropDbl d3alphar_dDelta2_dTau (void)
 Return the term \( \alpha^r_{\delta\delta\tau} \).
 
CoolPropDbl d3alphar_dDelta_dTau2 (void)
 Return the term \( \alpha^r_{\delta\tau\tau} \).
 
CoolPropDbl d3alphar_dTau3 (void)
 Return the term \( \alpha^r_{\tau\tau\tau} \).
 
CoolPropDbl d4alphar_dDelta4 (void)
 Return the term \( \alpha^r_{\delta\delta\delta\delta} \).
 
CoolPropDbl d4alphar_dDelta3_dTau (void)
 Return the term \( \alpha^r_{\delta\delta\delta\tau} \).
 
CoolPropDbl d4alphar_dDelta2_dTau2 (void)
 Return the term \( \alpha^r_{\delta\delta\tau\tau} \).
 
CoolPropDbl d4alphar_dDelta_dTau3 (void)
 Return the term \( \alpha^r_{\delta\tau\tau\tau} \).
 
CoolPropDbl d4alphar_dTau4 (void)
 Return the term \( \alpha^r_{\tau\tau\tau\tau} \).
 

Static Public Member Functions

static AbstractStatefactory (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. More...
 
static AbstractStatefactory (const std::string &backend, const std::vector< std::string > &fluid_names)
 A factory function to return a pointer to a new-allocated instance of one of the backends. More...
 

Protected Member Functions

bool isSupercriticalPhase (void)
 
bool isHomogeneousPhase (void)
 
bool isTwoPhase (void)
 
virtual CoolPropDbl calc_hmolar (void)
 Using this backend, calculate the molar enthalpy in J/mol.
 
virtual CoolPropDbl calc_smolar (void)
 Using this backend, calculate the molar entropy in J/mol/K.
 
virtual CoolPropDbl calc_umolar (void)
 Using this backend, calculate the molar internal energy in J/mol.
 
virtual CoolPropDbl calc_cpmolar (void)
 Using this backend, calculate the molar constant-pressure specific heat in J/mol/K.
 
virtual CoolPropDbl calc_cpmolar_idealgas (void)
 Using this backend, calculate the ideal gas molar constant-pressure specific heat in J/mol/K.
 
virtual CoolPropDbl calc_cvmolar (void)
 Using this backend, calculate the molar constant-volume specific heat in J/mol/K.
 
virtual CoolPropDbl calc_gibbsmolar (void)
 Using this backend, calculate the molar Gibbs function in J/mol.
 
virtual CoolPropDbl calc_helmholtzmolar (void)
 Using this backend, calculate the molar Helmholtz energy in J/mol.
 
virtual CoolPropDbl calc_speed_sound (void)
 Using this backend, calculate the speed of sound in m/s.
 
virtual CoolPropDbl calc_isothermal_compressibility (void)
 Using this backend, calculate the isothermal compressibility \( \kappa = -\frac{1}{v}\left.\frac{\partial v}{\partial p}\right|_T=\frac{1}{\rho}\left.\frac{\partial \rho}{\partial p}\right|_T\) in 1/Pa.
 
virtual CoolPropDbl calc_isobaric_expansion_coefficient (void)
 Using this backend, calculate the isobaric expansion coefficient \( \beta = \frac{1}{v}\left.\frac{\partial v}{\partial T}\right|_p = -\frac{1}{\rho}\left.\frac{\partial \rho}{\partial T}\right|_p\) in 1/K.
 
virtual CoolPropDbl calc_isentropic_expansion_coefficient (void)
 Using this backend, calculate the isentropic expansion coefficient \( \kappa_s = -\frac{c_p}{c_v}\frac{v}{p}\left.\frac{\partial p}{\partial v}\right|_T = \frac{\rho}{p}\left.\frac{\partial p}{\partial \rho}\right|_s\).
 
virtual CoolPropDbl calc_viscosity (void)
 Using this backend, calculate the viscosity in Pa-s.
 
virtual CoolPropDbl calc_conductivity (void)
 Using this backend, calculate the thermal conductivity in W/m/K.
 
virtual CoolPropDbl calc_surface_tension (void)
 Using this backend, calculate the surface tension in N/m.
 
virtual CoolPropDbl calc_molar_mass (void)
 Using this backend, calculate the molar mass in kg/mol.
 
virtual CoolPropDbl calc_acentric_factor (void)
 Using this backend, calculate the acentric factor.
 
virtual CoolPropDbl calc_pressure (void)
 Using this backend, calculate the pressure in Pa.
 
virtual CoolPropDbl calc_gas_constant (void)
 Using this backend, calculate the universal gas constant \(R_u\) in J/mol/K.
 
virtual CoolPropDbl calc_fugacity_coefficient (std::size_t i)
 Using this backend, calculate the fugacity coefficient (dimensionless)
 
virtual CoolPropDbl calc_fugacity (std::size_t i)
 Using this backend, calculate the fugacity in Pa.
 
virtual CoolPropDbl calc_chemical_potential (std::size_t i)
 Using this backend, calculate the chemical potential in J/mol.
 
virtual CoolPropDbl calc_PIP (void)
 Using this backend, calculate the phase identification parameter (PIP)
 
virtual void calc_excess_properties (void)
 Using this backend, calculate and cache the excess properties.
 
virtual CoolPropDbl calc_alphar (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r\) (dimensionless)
 
virtual CoolPropDbl calc_dalphar_dDelta (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta}\) (dimensionless)
 
virtual CoolPropDbl calc_dalphar_dTau (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d2alphar_dDelta2 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\delta}\) (dimensionless)
 
virtual CoolPropDbl calc_d2alphar_dDelta_dTau (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d2alphar_dTau2 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\tau\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d3alphar_dDelta3 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\delta\delta}\) (dimensionless)
 
virtual CoolPropDbl calc_d3alphar_dDelta2_dTau (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\delta\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d3alphar_dDelta_dTau2 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\tau\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d3alphar_dTau3 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\tau\tau\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d4alphar_dDelta4 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\delta\delta\delta}\) (dimensionless)
 
virtual CoolPropDbl calc_d4alphar_dDelta3_dTau (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\delta\delta\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d4alphar_dDelta2_dTau2 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\delta\tau\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d4alphar_dDelta_dTau3 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\delta\tau\tau\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d4alphar_dTau4 (void)
 Using this backend, calculate the residual Helmholtz energy term \(\alpha^r_{\tau\tau\tau\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_alpha0 (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0\) (dimensionless)
 
virtual CoolPropDbl calc_dalpha0_dDelta (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\delta}\) (dimensionless)
 
virtual CoolPropDbl calc_dalpha0_dTau (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d2alpha0_dDelta_dTau (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\delta\delta}\) (dimensionless)
 
virtual CoolPropDbl calc_d2alpha0_dDelta2 (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\delta\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d2alpha0_dTau2 (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\tau\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d3alpha0_dDelta3 (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\delta\delta\delta}\) (dimensionless)
 
virtual CoolPropDbl calc_d3alpha0_dDelta2_dTau (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\delta\delta\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d3alpha0_dDelta_dTau2 (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\delta\tau\tau}\) (dimensionless)
 
virtual CoolPropDbl calc_d3alpha0_dTau3 (void)
 Using this backend, calculate the ideal-gas Helmholtz energy term \(\alpha^0_{\tau\tau\tau}\) (dimensionless)
 
virtual void calc_reducing_state (void)
 
virtual CoolPropDbl calc_Tmax (void)
 Using this backend, calculate the maximum temperature in K.
 
virtual CoolPropDbl calc_Tmin (void)
 Using this backend, calculate the minimum temperature in K.
 
virtual CoolPropDbl calc_pmax (void)
 Using this backend, calculate the maximum pressure in Pa.
 
virtual CoolPropDbl calc_GWP20 (void)
 Using this backend, calculate the 20-year global warming potential (GWP)
 
virtual CoolPropDbl calc_GWP100 (void)
 Using this backend, calculate the 100-year global warming potential (GWP)
 
virtual CoolPropDbl calc_GWP500 (void)
 Using this backend, calculate the 500-year global warming potential (GWP)
 
virtual CoolPropDbl calc_ODP (void)
 Using this backend, calculate the ozone depletion potential (ODP)
 
virtual CoolPropDbl calc_flame_hazard (void)
 Using this backend, calculate the flame hazard.
 
virtual CoolPropDbl calc_health_hazard (void)
 Using this backend, calculate the health hazard.
 
virtual CoolPropDbl calc_physical_hazard (void)
 Using this backend, calculate the physical hazard.
 
virtual CoolPropDbl calc_dipole_moment (void)
 Using this backend, calculate the dipole moment in C-m (1 D = 3.33564e-30 C-m)
 
virtual CoolPropDbl calc_first_partial_deriv (parameters Of, parameters Wrt, parameters Constant)
 Calculate the first partial derivative for the desired derivative.
 
virtual CoolPropDbl calc_second_partial_deriv (parameters Of1, parameters Wrt1, parameters Constant1, parameters Wrt2, parameters Constant2)
 Calculate the second partial derivative using the given backend.
 
virtual CoolPropDbl calc_reduced_density (void)
 Using this backend, calculate the reduced density (rho/rhoc)
 
virtual CoolPropDbl calc_reciprocal_reduced_temperature (void)
 Using this backend, calculate the reciprocal reduced temperature (Tc/T)
 
virtual CoolPropDbl calc_Bvirial (void)
 Using this backend, calculate the second virial coefficient.
 
virtual CoolPropDbl calc_Cvirial (void)
 Using this backend, calculate the third virial coefficient.
 
virtual CoolPropDbl calc_dBvirial_dT (void)
 Using this backend, calculate the derivative dB/dT.
 
virtual CoolPropDbl calc_dCvirial_dT (void)
 Using this backend, calculate the derivative dC/dT.
 
virtual CoolPropDbl calc_compressibility_factor (void)
 Using this backend, calculate the compressibility factor Z \( Z = p/(\rho R T) \).
 
virtual std::string calc_name (void)
 Using this backend, get the name of the fluid.
 
virtual std::string calc_description (void)
 Using this backend, get the description of the fluid.
 
virtual CoolPropDbl calc_Ttriple (void)
 Using this backend, get the triple point temperature in K.
 
virtual CoolPropDbl calc_p_triple (void)
 Using this backend, get the triple point pressure in Pa.
 
virtual CoolPropDbl calc_T_critical (void)
 Using this backend, get the critical point temperature in K.
 
virtual CoolPropDbl calc_T_reducing (void)
 Using this backend, get the reducing point temperature in K.
 
virtual CoolPropDbl calc_p_critical (void)
 Using this backend, get the critical point pressure in Pa.
 
virtual CoolPropDbl calc_p_reducing (void)
 Using this backend, get the reducing point pressure in Pa.
 
virtual CoolPropDbl calc_rhomolar_critical (void)
 Using this backend, get the critical point molar density in mol/m^3.
 
virtual CoolPropDbl calc_rhomass_critical (void)
 Using this backend, get the critical point mass density in kg/m^3 - Added for IF97Backend which is mass based.
 
virtual CoolPropDbl calc_rhomolar_reducing (void)
 Using this backend, get the reducing point molar density in mol/m^3.
 
virtual void calc_phase_envelope (const std::string &type)
 Using this backend, construct the phase envelope, the variable type describes the type of phase envelope to be built.
 
virtual CoolPropDbl calc_rhomass (void)
 
virtual CoolPropDbl calc_hmass (void)
 
virtual CoolPropDbl calc_hmass_excess (void)
 
virtual CoolPropDbl calc_smass (void)
 
virtual CoolPropDbl calc_smass_excess (void)
 
virtual CoolPropDbl calc_cpmass (void)
 
virtual CoolPropDbl calc_cp0mass (void)
 
virtual CoolPropDbl calc_cvmass (void)
 
virtual CoolPropDbl calc_umass (void)
 
virtual CoolPropDbl calc_umass_excess (void)
 
virtual CoolPropDbl calc_gibbsmass (void)
 
virtual CoolPropDbl calc_gibbsmass_excess (void)
 
virtual CoolPropDbl calc_helmholtzmass (void)
 
virtual CoolPropDbl calc_helmholtzmass_excess (void)
 
virtual CoolPropDbl calc_volumemass_excess (void)
 
virtual void update_states (void)
 Update the states after having changed the reference state for enthalpy and entropy.
 
virtual CoolPropDbl calc_melting_line (int param, int given, CoolPropDbl value)
 
virtual CoolPropDbl calc_saturation_ancillary (parameters param, int Q, parameters given, double value)
 
virtual phases calc_phase (void)
 Using this backend, calculate the phase.
 
virtual void calc_specify_phase (phases phase)
 Using this backend, specify the phase to be used for all further calculations.
 
virtual void calc_unspecify_phase (void)
 Using this backend, unspecify the phase.
 
virtual std::vector< std::string > calc_fluid_names (void)
 Using this backend, get a vector of fluid names.
 
virtual const CoolProp::SimpleStatecalc_state (const std::string &state)
 Using this backend, calculate a phase given by the state string. More...
 
virtual const CoolProp::PhaseEnvelopeDatacalc_phase_envelope_data (void)
 
virtual std::vector< CoolPropDbl > calc_mole_fractions_liquid (void)
 
virtual std::vector< CoolPropDbl > calc_mole_fractions_vapor (void)
 
virtual const std::vector< CoolPropDbl > calc_mass_fractions (void)
 
virtual CoolPropDbl calc_fraction_min (void)
 Get the minimum fraction (mole, mass, volume) for incompressible fluid.
 
virtual CoolPropDbl calc_fraction_max (void)
 Get the maximum fraction (mole, mass, volume) for incompressible fluid.
 
virtual CoolPropDbl calc_T_freeze (void)
 
virtual CoolPropDbl calc_first_saturation_deriv (parameters Of1, parameters Wrt1)
 
virtual CoolPropDbl calc_second_saturation_deriv (parameters Of1, parameters Wrt1, parameters Wrt2)
 
virtual CoolPropDbl calc_first_two_phase_deriv (parameters Of, parameters Wrt, parameters Constant)
 
virtual CoolPropDbl calc_second_two_phase_deriv (parameters Of, parameters Wrt, parameters Constant, parameters Wrt2, parameters Constant2)
 
virtual CoolPropDbl calc_first_two_phase_deriv_splined (parameters Of, parameters Wrt, parameters Constant, CoolPropDbl x_end)
 
virtual CoolPropDbl calc_saturated_liquid_keyed_output (parameters key)
 
virtual CoolPropDbl calc_saturated_vapor_keyed_output (parameters key)
 
virtual void calc_ideal_curve (const std::string &type, std::vector< double > &T, std::vector< double > &p)
 
virtual CoolPropDbl calc_T (void)
 Using this backend, get the temperature.
 
virtual CoolPropDbl calc_rhomolar (void)
 Using this backend, get the molar density in mol/m^3.
 
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.
 
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.
 
virtual void calc_conformal_state (const std::string &reference_fluid, CoolPropDbl &T, CoolPropDbl &rhomolar)
 
virtual void calc_viscosity_contributions (CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical)
 
virtual void calc_conductivity_contributions (CoolPropDbl &dilute, CoolPropDbl &initial_density, CoolPropDbl &residual, CoolPropDbl &critical)
 
virtual std::vector< CriticalStatecalc_all_critical_points (void)
 
virtual void calc_build_spinodal ()
 
virtual SpinodalData calc_get_spinodal_data ()
 
virtual void calc_criticality_contour_values (double &L1star, double &M1star)
 
virtual void mass_to_molar_inputs (CoolProp::input_pairs &input_pair, CoolPropDbl &value1, CoolPropDbl &value2)
 Convert mass-based input pair to molar-based input pair; If molar-based, do nothing. More...
 
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.
 

Protected Attributes

long _fluid_type
 Some administrative variables.
 
phases _phase
 The key for the phase from CoolProp::phases enum.
 
phases imposed_phase_index
 If the phase is imposed, the imposed phase index.
 
SimpleState _critical
 Two important points.
 
SimpleState _reducing
 
CachedElement _molar_mass
 Molar mass [mol/kg].
 
CachedElement _gas_constant
 Universal gas constant [J/mol/K].
 
double _rhomolar
 Bulk values.
 
double _T
 
double _p
 
double _Q
 
double _R
 
CachedElement _tau
 
CachedElement _delta
 
CachedElement _viscosity
 Transport properties.
 
CachedElement _conductivity
 
CachedElement _surface_tension
 
CachedElement _hmolar
 
CachedElement _smolar
 
CachedElement _umolar
 
CachedElement _logp
 
CachedElement _logrhomolar
 
CachedElement _cpmolar
 
CachedElement _cp0molar
 
CachedElement _cvmolar
 
CachedElement _speed_sound
 
CachedElement _gibbsmolar
 
CachedElement _helmholtzmolar
 
CachedElement _hmolar_excess
 Excess properties.
 
CachedElement _smolar_excess
 
CachedElement _gibbsmolar_excess
 
CachedElement _umolar_excess
 
CachedElement _volumemolar_excess
 
CachedElement _helmholtzmolar_excess
 
CachedElement _rhoLanc
 Ancillary values.
 
CachedElement _rhoVanc
 
CachedElement _pLanc
 
CachedElement _pVanc
 
CachedElement _TLanc
 
CachedElement _TVanc
 
CachedElement _fugacity_coefficient
 
CachedElement _rho_spline
 Smoothing values.
 
CachedElement _drho_spline_dh__constp
 
CachedElement _drho_spline_dp__consth
 
CachedElement _alpha0
 Cached low-level elements for in-place calculation of other properties.
 
CachedElement _dalpha0_dTau
 
CachedElement _dalpha0_dDelta
 
CachedElement _d2alpha0_dTau2
 
CachedElement _d2alpha0_dDelta_dTau
 
CachedElement _d2alpha0_dDelta2
 
CachedElement _d3alpha0_dTau3
 
CachedElement _d3alpha0_dDelta_dTau2
 
CachedElement _d3alpha0_dDelta2_dTau
 
CachedElement _d3alpha0_dDelta3
 
CachedElement _alphar
 
CachedElement _dalphar_dTau
 
CachedElement _dalphar_dDelta
 
CachedElement _d2alphar_dTau2
 
CachedElement _d2alphar_dDelta_dTau
 
CachedElement _d2alphar_dDelta2
 
CachedElement _d3alphar_dTau3
 
CachedElement _d3alphar_dDelta_dTau2
 
CachedElement _d3alphar_dDelta2_dTau
 
CachedElement _d3alphar_dDelta3
 
CachedElement _d4alphar_dTau4
 
CachedElement _d4alphar_dDelta_dTau3
 
CachedElement _d4alphar_dDelta2_dTau2
 
CachedElement _d4alphar_dDelta3_dTau
 
CachedElement _d4alphar_dDelta4
 
CachedElement _dalphar_dDelta_lim
 
CachedElement _d2alphar_dDelta2_lim
 
CachedElement _d2alphar_dDelta_dTau_lim
 
CachedElement _d3alphar_dDelta2_dTau_lim
 
CachedElement _rhoLmolar
 Two-Phase variables.
 
CachedElement _rhoVmolar
 

Member Function Documentation

§ backend_name()

virtual std::string CoolProp::AbstractState::backend_name ( void  )
pure virtual

Get a string representation of the backend - for instance "HelmholtzEOSMixtureBackend" for the core mixture model in CoolProp.

Must be overloaded by the backend to provide the backend's name

Implemented in CoolProp::PengRobinsonBackend, CoolProp::SRKBackend, CoolProp::BicubicBackend, CoolProp::HelmholtzEOSMixtureBackend, CoolProp::VTPRBackend, CoolProp::HelmholtzEOSBackend, CoolProp::REFPROPMixtureBackend, CoolProp::IncompressibleBackend, CoolProp::REFPROPBackend, CoolProp::IF97Backend, and CoolProp::TTSEBackend.

§ build_phase_envelope()

void CoolProp::AbstractState::build_phase_envelope ( const std::string &  type = "")

Construct the phase envelope for a mixture.

Parameters
typecurrently a dummy variable that is not used

§ calc_saturation_ancillary()

virtual CoolPropDbl CoolProp::AbstractState::calc_saturation_ancillary ( parameters  param,
int  Q,
parameters  given,
double  value 
)
inlineprotectedvirtual
Parameters
paramThe key for the parameter to be returned
QThe quality for the parameter that is given (0 = saturated liquid, 1 = saturated vapor)
givenThe key for the parameter that is given
valueThe value for the parameter that is given

Reimplemented in CoolProp::HelmholtzEOSMixtureBackend.

§ calc_state()

virtual const CoolProp::SimpleState& CoolProp::AbstractState::calc_state ( const std::string &  state)
inlineprotectedvirtual

Using this backend, calculate a phase given by the state string.

Parameters
stateA string that describes the state desired, one of "hs_anchor", "critical"/"crit", "reducing"

Reimplemented in CoolProp::HelmholtzEOSMixtureBackend.

§ change_EOS()

void CoolProp::AbstractState::change_EOS ( const std::size_t  i,
const std::string &  EOS_name 
)
inline

Change the equation of state for a given component to a specified EOS.

Parameters
iIndex of the component to change (if a pure fluid, i=0)
EOS_nameName of the EOS to use (something like "SRK", "PR", "XiangDeiters", but backend-specific)
Note
Calls the calc_change_EOS function of the implementation

§ clear()

bool CoolProp::AbstractState::clear ( )
virtual

Clear all the cached values.

Ancillary curve values

Bulk values

Smoothing values

Cached low-level elements for in-place calculation of other properties

Two-Phase variables

Transport properties

Reimplemented in CoolProp::HelmholtzEOSMixtureBackend, CoolProp::IncompressibleBackend, and CoolProp::IF97Backend.

§ conformal_state()

void CoolProp::AbstractState::conformal_state ( const std::string &  reference_fluid,
CoolPropDbl &  T,
CoolPropDbl &  rhomolar 
)
inline

Find the conformal state needed for ECS.

Parameters
reference_fluidThe reference fluid for which the conformal state will be calculated
TTemperature (initial guess must be provided, or < 0 to start with unity shape factors)
rhomolarMolar density (initial guess must be provided, or < 0 to start with unity shape factors)

§ factory() [1/2]

static AbstractState* CoolProp::AbstractState::factory ( const std::string &  backend,
const std::string &  fluid_names 
)
inlinestatic

A factory function to return a pointer to a new-allocated instance of one of the backends.

This is a convenience function to allow for the use of '&' delimited fluid names. Slightly less computationally efficient than the

Parameters
backendThe backend in use, one of "HEOS", "REFPROP", etc.
fluid_namesFluid names as a '&' delimited string
Returns

§ factory() [2/2]

AbstractState * CoolProp::AbstractState::factory ( const std::string &  backend,
const std::vector< std::string > &  fluid_names 
)
static

A factory function to return a pointer to a new-allocated instance of one of the backends.

Parameters
backendThe backend in use, "HEOS", "REFPROP", etc.
fluid_namesA vector of strings of the fluid names
Returns
A pointer to the instance generated

Several backends are possible:

  1. "?" : The backend is unknown, we will parse the fluid string to determine the backend to be used. Probably will use HEOS backend (see below)
  2. "HEOS" : The Helmholtz Equation of State backend for use with pure and pseudo-pure fluids, and mixtures, all of which are based on multi-parameter Helmholtz Energy equations of state. The fluid part of the string should then either be
    1. A pure or pseudo-pure fluid name (eg. "PROPANE" or "R410A"), yielding a HelmholtzEOSBackend instance.
    2. A string that encodes the components of the mixture with a "&" between them (e.g. "R32&R125"), yielding a HelmholtzEOSMixtureBackend instance.
  3. "REFPROP" : The REFPROP backend will be used. The fluid part of the string should then either be
    1. A pure or pseudo-pure fluid name (eg. "PROPANE" or "R410A"), yielding a REFPROPBackend instance.
    2. A string that encodes the components of the mixture with a "&" between them (e.g. "R32&R125"), yielding a REFPROPMixtureBackend instance.
  4. "INCOMP": The incompressible backend will be used
  5. "TTSE&XXXX": The TTSE backend will be used, and the tables will be generated using the XXXX backend where XXXX is one of the base backends("HEOS", "REFPROP", etc. )
  6. "BICUBIC&XXXX": The Bicubic backend will be used, and the tables will be generated using the XXXX backend where XXXX is one of the base backends("HEOS", "REFPROP", etc. )

Very Important!! : Use a smart pointer to manage the pointer returned. In older versions of C++, you can use std::tr1::smart_ptr. In C++2011 you can use std::shared_ptr

§ first_partial_deriv()

CoolPropDbl CoolProp::AbstractState::first_partial_deriv ( parameters  Of,
parameters  Wrt,
parameters  Constant 
)
inline

The first partial derivative in homogeneous phases.

\[ \left(\frac{\partial A}{\partial B}\right)_C = \frac{\left(\frac{\partial A}{\partial \tau}\right)_\delta\left(\frac{\partial C}{\partial \delta}\right)_\tau-\left(\frac{\partial A}{\partial \delta}\right)_\tau\left(\frac{\partial C}{\partial \tau}\right)_\delta}{\left(\frac{\partial B}{\partial \tau}\right)_\delta\left(\frac{\partial C}{\partial \delta}\right)_\tau-\left(\frac{\partial B}{\partial \delta}\right)_\tau\left(\frac{\partial C}{\partial \tau}\right)_\delta} = \frac{N}{D}\]

§ first_saturation_deriv()

CoolPropDbl CoolProp::AbstractState::first_saturation_deriv ( parameters  Of1,
parameters  Wrt1 
)
inline

The first partial derivative along the saturation curve.

Implementing the algorithms and ideas of: Matthis Thorade, Ali Saadat, "Partial derivatives of thermodynamic state properties for dynamic simulation", Environmental Earth Sciences, December 2013, Volume 70, Issue 8, pp 3497-3503

Basically the idea is that the p-T derivative is given by Clapeyron relations:

\[ \left(\frac{\partial T}{\partial p}\right)_{\sigma} = T\left(\frac{v'' - v'}{h'' - h'}\right)_{\sigma} \]

and then other derivatives can be obtained along the saturation curve from

\[ \left(\frac{\partial y}{\partial p}\right)_{\sigma} = \left(\frac{\partial y}{\partial p}\right)+\left(\frac{\partial y}{\partial T}\right)\left(\frac{\partial T}{\partial p}\right)_{\sigma} \]

\[ \left(\frac{\partial y}{\partial T}\right)_{\sigma} = \left(\frac{\partial y}{\partial T}\right)+\left(\frac{\partial y}{\partial p}\right)\left(\frac{\partial p}{\partial T}\right)_{\sigma} \]

where derivatives without the \( \sigma \) are homogeneous (conventional) derivatives.

Parameters
Of1The parameter that the derivative is taken of
Wrt1The parameter that the derivative is taken with respect to

§ first_two_phase_deriv()

double CoolProp::AbstractState::first_two_phase_deriv ( parameters  Of,
parameters  Wrt,
parameters  Constant 
)
inline

Calculate the first "two-phase" derivative as described by Thorade and Sadaat, EAS, 2013.

Implementing the algorithms and ideas of: Matthis Thorade, Ali Saadat, "Partial derivatives of thermodynamic state properties for dynamic simulation", Environmental Earth Sciences, December 2013, Volume 70, Issue 8, pp 3497-3503

Spline evaluation is as described in: S Quoilin, I Bell, A Desideri, P Dewallef, V Lemort, "Methods to increase the robustness of finite-volume flow models in thermodynamic systems", Energies 7 (3), 1621-1640

Note
Not all derivatives are supported!
Parameters
OfThe parameter to be derived
WrtThe parameter that the derivative is taken with respect to
ConstantThe parameter that is held constant
Returns

§ first_two_phase_deriv_splined()

double CoolProp::AbstractState::first_two_phase_deriv_splined ( parameters  Of,
parameters  Wrt,
parameters  Constant,
double  x_end 
)
inline

Calculate the first "two-phase" derivative as described by Thorade and Sadaat, EAS, 2013.

Implementing the algorithms and ideas of: Matthis Thorade, Ali Saadat, "Partial derivatives of thermodynamic state properties for dynamic simulation", Environmental Earth Sciences, December 2013, Volume 70, Issue 8, pp 3497-3503

Spline evaluation is as described in: S Quoilin, I Bell, A Desideri, P Dewallef, V Lemort, "Methods to increase the robustness of finite-volume flow models in thermodynamic systems", Energies 7 (3), 1621-1640

Note
Not all derivatives are supported! If you need all three currently supported values (drho_dh__p, drho_dp__h and rho_spline), you should calculate drho_dp__h first to avoid duplicate calculations.
Parameters
OfThe parameter to be derived
WrtThe parameter that the derivative is taken with respect to
ConstantThe parameter that is held constant
x_endThe end vapor quality at which the spline is defined (spline is active in [0, x_end])
Returns

§ fundamental_derivative_of_gas_dynamics()

double CoolProp::AbstractState::fundamental_derivative_of_gas_dynamics ( void  )

Return the fundamental derivative of gas dynamics \( \Gamma \).

see also Colonna et al, FPE, 2010

\[ \Gamma = 1+\frac{\rho}{c}\left(\frac{partial c}{\partial \rho}\right)_{s} = 1+\frac{\rho}{2c^2}\left(\frac{partial^2 p}{\partial \rho^2}\right)_{s} = 1+\frac{v^3}{2c^2}\left(\frac{partial^2 p}{\partial v^2}\right)_{s}\]

Note: densities are mass-based densities, not mole-based densities

§ get_fluid_constant()

virtual const double CoolProp::AbstractState::get_fluid_constant ( std::size_t  i,
parameters  param 
) const
inlinevirtual

Get a constant for one of the fluids forming this mixture.

Parameters
iIndex (0-based) of the fluid
paramparameter you want to obtain (probably one that is a trivial parameter)

Reimplemented in CoolProp::HelmholtzEOSMixtureBackend, and CoolProp::AbstractCubicBackend.

§ get_reducing_state()

virtual const CoolProp::SimpleState& CoolProp::AbstractState::get_reducing_state ( )
inlinevirtual

Get the state that is used in the equation of state or mixture model to reduce the state.

For pure fluids this is usually, but not always, the critical point. For mixture models, it is usually composition dependent

Reimplemented in CoolProp::HelmholtzEOSMixtureBackend.

§ ideal_curve()

void CoolProp::AbstractState::ideal_curve ( const std::string &  type,
std::vector< double > &  T,
std::vector< double > &  p 
)
inline

Calculate an ideal curve for a pure fluid.

Parameters
typeThe type of ideal curve you would like to calculate - "Ideal", "Boyle", "Joule-Thomson", "Joule Inversion", etc.
TThe temperatures along the curve in K
pThe pressures along the curve in Pa

§ mass_to_molar_inputs()

void CoolProp::AbstractState::mass_to_molar_inputs ( CoolProp::input_pairs input_pair,
CoolPropDbl &  value1,
CoolPropDbl &  value2 
)
protectedvirtual

Convert mass-based input pair to molar-based input pair; If molar-based, do nothing.

< Mass density in kg/m^3, Temperature in K

< Entropy in J/kg/K, Temperature in K

< Mass density in kg/m^3, Pressure in Pa

< Mass density in kg/m^3, molar quality

< Enthalpy in J/kg, Pressure in Pa

< Pressure in Pa, Entropy in J/kg/K

< Pressure in Pa, Internal energy in J/kg

< Enthalpy in J/kg, Entropy in J/kg/K

< Entropy in J/kg/K, Internal energy in J/kg

< Mass density in kg/m^3, Enthalpy in J/kg

< Mass density in kg/m^3, Entropy in J/kg/K

< Mass density in kg/m^3, Internal energy in J/kg

§ melting_line()

double CoolProp::AbstractState::melting_line ( int  param,
int  given,
double  value 
)

Return a value from the melting line.

Parameters
paramThe key for the parameter to be returned
givenThe key for the parameter that is given
valueThe value for the parameter that is given

§ saturation_ancillary()

double CoolProp::AbstractState::saturation_ancillary ( parameters  param,
int  Q,
parameters  given,
double  value 
)

Return the value from a saturation ancillary curve (if the backend implements it)

Parameters
paramThe key for the parameter to be returned
QThe quality for the parameter that is given (0 = saturated liquid, 1 = saturated vapor)
givenThe key for the parameter that is given
valueThe value for the parameter that is given

§ second_partial_deriv()

CoolPropDbl CoolProp::AbstractState::second_partial_deriv ( parameters  Of1,
parameters  Wrt1,
parameters  Constant1,
parameters  Wrt2,
parameters  Constant2 
)
inline

The second partial derivative in homogeneous phases.

The first partial derivative (CoolProp::AbstractState::first_partial_deriv) can be expressed as

\[ \left(\frac{\partial A}{\partial B}\right)_C = \frac{\left(\frac{\partial A}{\partial T}\right)_\rho\left(\frac{\partial C}{\partial \rho}\right)_T-\left(\frac{\partial A}{\partial \rho}\right)_T\left(\frac{\partial C}{\partial T}\right)_\rho}{\left(\frac{\partial B}{\partial T}\right)_\rho\left(\frac{\partial C}{\partial \rho}\right)_T-\left(\frac{\partial B}{\partial \rho}\right)_T\left(\frac{\partial C}{\partial T}\right)_\rho} = \frac{N}{D}\]

and the second derivative can be expressed as

\[ \frac{\partial}{\partial D}\left(\left(\frac{\partial A}{\partial B}\right)_C\right)_E = \frac{\frac{\partial}{\partial T}\left( \left(\frac{\partial A}{\partial B}\right)_C \right)_\rho\left(\frac{\partial E}{\partial \rho}\right)_T-\frac{\partial}{\partial \rho}\left(\left(\frac{\partial A}{\partial B}\right)_C\right)_T\left(\frac{\partial E}{\partial T}\right)_\rho}{\left(\frac{\partial D}{\partial T}\right)_\rho\left(\frac{\partial E}{\partial \rho}\right)_T-\left(\frac{\partial D}{\partial \rho}\right)_T\left(\frac{\partial E}{\partial T}\right)_\rho} \]

which can be expressed in parts as

\[\left(\frac{\partial N}{\partial \rho}\right)_{T} = \left(\frac{\partial A}{\partial T}\right)_\rho\left(\frac{\partial^2 C}{\partial \rho^2}\right)_{T}+\left(\frac{\partial^2 A}{\partial T\partial\rho}\right)\left(\frac{\partial C}{\partial \rho}\right)_{T}-\left(\frac{\partial A}{\partial \rho}\right)_T\left(\frac{\partial^2 C}{\partial T\partial\rho}\right)-\left(\frac{\partial^2 A}{\partial \rho^2}\right)_{T}\left(\frac{\partial C}{\partial T}\right)_\rho\]

\[\left(\frac{\partial D}{\partial \rho}\right)_{T} = \left(\frac{\partial B}{\partial T}\right)_\rho\left(\frac{\partial^2 C}{\partial \rho^2}\right)_{T}+\left(\frac{\partial^2 B}{\partial T\partial\rho}\right)\left(\frac{\partial C}{\partial \rho}\right)_{T}-\left(\frac{\partial B}{\partial \rho}\right)_T\left(\frac{\partial^2 C}{\partial T\partial\rho}\right)-\left(\frac{\partial^2 B}{\partial \rho^2}\right)_{T}\left(\frac{\partial C}{\partial T}\right)_\rho\]

\[\left(\frac{\partial N}{\partial T}\right)_{\rho} = \left(\frac{\partial A}{\partial T}\right)_\rho\left(\frac{\partial^2 C}{\partial \rho\partial T}\right)+\left(\frac{\partial^2 A}{\partial T^2}\right)_\rho\left(\frac{\partial C}{\partial \rho}\right)_{T}-\left(\frac{\partial A}{\partial \rho}\right)_T\left(\frac{\partial^2 C}{\partial T^2}\right)_\rho-\left(\frac{\partial^2 A}{\partial \rho\partial T}\right)\left(\frac{\partial C}{\partial T}\right)_\rho\]

\[\left(\frac{\partial D}{\partial T}\right)_{\rho} = \left(\frac{\partial B}{\partial T}\right)_\rho\left(\frac{\partial^2 C}{\partial \rho\partial T}\right)+\left(\frac{\partial^2 B}{\partial T^2}\right)_\rho\left(\frac{\partial C}{\partial \rho}\right)_{T}-\left(\frac{\partial B}{\partial \rho}\right)_T\left(\frac{\partial^2 C}{\partial T^2}\right)_\rho-\left(\frac{\partial^2 B}{\partial \rho\partial T}\right)\left(\frac{\partial C}{\partial T}\right)_\rho\]

\[\frac{\partial}{\partial \rho}\left( \left(\frac{\partial A}{\partial B}\right)_C \right)_T = \frac{D\left(\frac{\partial N}{\partial \rho}\right)_{T}-N\left(\frac{\partial D}{\partial \rho}\right)_{\tau}}{D^2}\]

\[\frac{\partial}{\partial T}\left( \left(\frac{\partial A}{\partial B}\right)_C \right)_\rho = \frac{D\left(\frac{\partial N}{\partial T}\right)_{\rho}-N\left(\frac{\partial D}{\partial T}\right)_{\rho}}{D^2}\]

The terms \( N \) and \( D \) are the numerator and denominator from CoolProp::AbstractState::first_partial_deriv respectively

§ second_saturation_deriv()

CoolPropDbl CoolProp::AbstractState::second_saturation_deriv ( parameters  Of1,
parameters  Wrt1,
parameters  Wrt2 
)
inline

The second partial derivative along the saturation curve.

Implementing the algorithms and ideas of: Matthis Thorade, Ali Saadat, "Partial derivatives of thermodynamic state properties for dynamic simulation", Environmental Earth Sciences, December 2013, Volume 70, Issue 8, pp 3497-3503

Like with first_saturation_deriv, we can express the derivative as

\[ \left(\frac{\partial y}{\partial T}\right)_{\sigma} = \left(\frac{\partial y}{\partial T}\right)+\left(\frac{\partial y}{\partial p}\right)\left(\frac{\partial p}{\partial T}\right)_{\sigma} \]

where \( y \) is already a saturation derivative. So you might end up with something like

\[ \left(\frac{\partial \left(\frac{\partial T}{\partial p}\right)_{\sigma}}{\partial T}\right)_{\sigma} = \left(\frac{\partial \left(\frac{\partial T}{\partial p}\right)_{\sigma}}{\partial T}\right)+\left(\frac{\partial \left(\frac{\partial T}{\partial p}\right)_{\sigma}}{\partial p}\right)\left(\frac{\partial p}{\partial T}\right)_{\sigma} \]

Parameters
Of1The parameter that the first derivative is taken of
Wrt1The parameter that the first derivative is taken with respect to
Wrt2The parameter that the second derivative is taken with respect to

§ second_two_phase_deriv()

double CoolProp::AbstractState::second_two_phase_deriv ( parameters  Of,
parameters  Wrt1,
parameters  Constant1,
parameters  Wrt2,
parameters  Constant2 
)
inline

Calculate the second "two-phase" derivative as described by Thorade and Sadaat, EAS, 2013.

Implementing the algorithms and ideas of: Matthis Thorade, Ali Saadat, "Partial derivatives of thermodynamic state properties for dynamic simulation", Environmental Earth Sciences, December 2013, Volume 70, Issue 8, pp 3497-3503

Note
Not all derivatives are supported!
Parameters
OfThe parameter to be derived
Wrt1The parameter that the derivative is taken with respect to in the first derivative
Constant1The parameter that is held constant in the first derivative
Wrt2The parameter that the derivative is taken with respect to in the second derivative
Constant2The parameter that is held constant in the second derivative
Returns

§ set_reference_stateD()

virtual void CoolProp::AbstractState::set_reference_stateD ( double  T,
double  rhomolar,
double  hmolar0,
double  smolar0 
)
inlinevirtual

Set the reference state based on a thermodynamic state point specified by temperature and molar density.

Parameters
TTemperature at reference state [K]
rhomolarMolar density at reference state [mol/m^3]
hmolar0Molar enthalpy at reference state [J/mol]
smolar0Molar entropy at reference state [J/mol/K]

Reimplemented in CoolProp::HelmholtzEOSMixtureBackend.

§ set_reference_stateS()

virtual void CoolProp::AbstractState::set_reference_stateS ( const std::string &  reference_state)
inlinevirtual

Set the reference state based on a string representation.

Parameters
reference_stateThe reference state to use, one of
Reference State Description
"IIR" h = 200 kJ/kg, s=1 kJ/kg/K at 0C saturated liquid
"ASHRAE" h = 0, s = 0 @ -40C saturated liquid
"NBP" h = 0, s = 0 @ 1.0 bar saturated liquid
"DEF" Reset to the default reference state for the fluid
"RESET" Remove the offset

The offset in the ideal gas Helmholtz energy can be obtained from

\[ \displaystyle\frac{\Delta s}{R_u/M}+\frac{\Delta h}{(R_u/M)T}\tau \]

where \( \Delta s = s-s_{spec} \) and \( \Delta h = h-h_{spec} \)

Reimplemented in CoolProp::HelmholtzEOSMixtureBackend.

§ tangent_plane_distance()

double CoolProp::AbstractState::tangent_plane_distance ( const double  T,
const double  p,
const std::vector< double > &  w,
const double  rhomolar_guess = -1 
)
inline

Return the tangent plane distance for a given trial composition w.

Parameters
TTemperature (K)
pPressure (Pa)
wThe trial composition
rhomolar_guess(mol/m^3) The molar density guess value (if <0 (default), not used; if >0, guess value will be used in flash evaluation)

\[ tpd(w) = \sum_i w_i(\ln w_i + \ln \phi_i(w) - d_i) \]

with

\[ d_i = \ln z_i + \ln \phi_i(z) \]

Or you can express the \( tpd \) in terms of fugacity (See Table 7.3 from GERG 2004 monograph) since \( \ln \phi_i = \ln f_i - \ln p -\ln z_i\) thus

\[ d_i = \ln f_i(z) - \ln p\]

and

\[ tpd(w) = \sum_i w_i(\ln f_i(w) - \ln p - d_i) \]

and the \( \ln p \) cancel, leaving

\[ tpd(w) = \sum_i w_i(\ln f_i(w) - \ln f_i(z)) \]


The documentation for this class was generated from the following files: