Skip to content

Commit

Permalink
Removed pEOS and replaced it with EOS() function - reference state wo…
Browse files Browse the repository at this point in the history
…rks properly; Closes #524
  • Loading branch information
ibell committed Mar 6, 2015
1 parent 7a9f98d commit 6b6f2db
Show file tree
Hide file tree
Showing 15 changed files with 375 additions and 270 deletions.
9 changes: 5 additions & 4 deletions include/CoolPropFluid.h
Expand Up @@ -45,7 +45,8 @@ struct CriticalRegionSplines{
std::vector<double> cL, cV;
bool enabled;
CriticalRegionSplines(){enabled = false;};
void get_densities(double T, double rho_min, double rho_crit, double rho_max, double &rhoL, double &rhoV){
const void get_densities(double T, double rho_min, double rho_crit, double rho_max, double &rhoL, double &rhoV) const
{
int Nsoln = -1, Ngood = 0;
double rho1 =0, rho2 = 0, rho3 = 0;

Expand Down Expand Up @@ -453,7 +454,7 @@ class CoolPropFluid {
public:
CoolPropFluid(){};
~CoolPropFluid(){};
EquationOfState *pEOS; ///< A pointer to the currently used EOS
EquationOfState &EOS() {return EOSVector[0];}; ///< Get a reference to the equation of state
std::vector<EquationOfState> EOSVector; ///< The equations of state that could be used for this fluid

std::string name; ///< The name of the fluid
Expand All @@ -469,8 +470,8 @@ class CoolPropFluid {
triple_liquid, ///< The saturated liquid state at the triple point temperature
triple_vapor; ///< The saturated vapor state at the triple point temperature

double gas_constant(){ return pEOS->R_u; };
double molar_mass(){ return pEOS->molar_mass; };
double gas_constant(){ return EOS().R_u; };
double molar_mass(){ return EOS().molar_mass; };
};


Expand Down
46 changes: 23 additions & 23 deletions src/Backends/Helmholtz/FlashRoutines.cpp
Expand Up @@ -149,8 +149,8 @@ void FlashRoutines::HQ_flash(HelmholtzEOSMixtureBackend &HEOS, CoolPropDbl Tgues
if (Tguess < 0){
options.use_guesses = true;
options.T = Tguess;
CoolProp::SaturationAncillaryFunction &rhoL = HEOS.get_components()[0]->ancillaries.rhoL;
CoolProp::SaturationAncillaryFunction &rhoV = HEOS.get_components()[0]->ancillaries.rhoV;
CoolProp::SaturationAncillaryFunction &rhoL = HEOS.get_components()[0].ancillaries.rhoL;
CoolProp::SaturationAncillaryFunction &rhoV = HEOS.get_components()[0].ancillaries.rhoV;
options.rhoL = rhoL.evaluate(Tguess);
options.rhoV = rhoV.evaluate(Tguess);
}
Expand Down Expand Up @@ -227,7 +227,7 @@ void FlashRoutines::QT_flash(HelmholtzEOSMixtureBackend &HEOS)
Tmin_sat = std::max(Tmin_satL, Tmin_satV) - 1e-13;

// Get a reference to keep the code a bit cleaner
CriticalRegionSplines &splines = HEOS.components[0]->pEOS->critical_region_splines;
const CriticalRegionSplines &splines = HEOS.components[0].EOS().critical_region_splines;

// If exactly(ish) at the critical temperature, liquid and vapor have the critial density
if ((get_config_bool(CRITICAL_WITHIN_1UK) && std::abs(T-Tmax_sat)< 1e-6) || std::abs(T-Tmax_sat)< 1e-12){
Expand All @@ -248,7 +248,7 @@ void FlashRoutines::QT_flash(HelmholtzEOSMixtureBackend &HEOS)
HEOS._p = 0.5*HEOS.SatV->p() + 0.5*HEOS.SatL->p();
HEOS._rhomolar = 1/(HEOS._Q/HEOS.SatV->rhomolar() + (1 - HEOS._Q)/HEOS.SatL->rhomolar());
}
else if (!(HEOS.components[0]->pEOS->pseudo_pure))
else if (!(HEOS.components[0].EOS().pseudo_pure))
{
// Set some imput options
SaturationSolvers::saturation_T_pure_Akasaka_options options;
Expand All @@ -263,11 +263,11 @@ void FlashRoutines::QT_flash(HelmholtzEOSMixtureBackend &HEOS)
else{
// Pseudo-pure fluid
CoolPropDbl rhoLanc = _HUGE, rhoVanc = _HUGE, rhoLsat = _HUGE, rhoVsat = _HUGE;
CoolPropDbl psatLanc = HEOS.components[0]->ancillaries.pL.evaluate(HEOS._T); // These ancillaries are used explicitly
CoolPropDbl psatVanc = HEOS.components[0]->ancillaries.pV.evaluate(HEOS._T); // These ancillaries are used explicitly
CoolPropDbl psatLanc = HEOS.components[0].ancillaries.pL.evaluate(HEOS._T); // These ancillaries are used explicitly
CoolPropDbl psatVanc = HEOS.components[0].ancillaries.pV.evaluate(HEOS._T); // These ancillaries are used explicitly
try{
rhoLanc = HEOS.components[0]->ancillaries.rhoL.evaluate(HEOS._T);
rhoVanc = HEOS.components[0]->ancillaries.rhoV.evaluate(HEOS._T);
rhoLanc = HEOS.components[0].ancillaries.rhoL.evaluate(HEOS._T);
rhoVanc = HEOS.components[0].ancillaries.rhoV.evaluate(HEOS._T);

if (!ValidNumber(rhoLanc) || !ValidNumber(rhoVanc))
{
Expand Down Expand Up @@ -329,14 +329,14 @@ void FlashRoutines::PQ_flash(HelmholtzEOSMixtureBackend &HEOS)
{
if (HEOS.is_pure_or_pseudopure)
{
if (HEOS.components[0]->pEOS->pseudo_pure){
if (HEOS.components[0].EOS().pseudo_pure){
// It is a pseudo-pure mixture

HEOS._TLanc = HEOS.components[0]->ancillaries.pL.invert(HEOS._p);
HEOS._TVanc = HEOS.components[0]->ancillaries.pV.invert(HEOS._p);
HEOS._TLanc = HEOS.components[0].ancillaries.pL.invert(HEOS._p);
HEOS._TVanc = HEOS.components[0].ancillaries.pV.invert(HEOS._p);
// Get guesses for the ancillaries for density
CoolPropDbl rhoL = HEOS.components[0]->ancillaries.rhoL.evaluate(HEOS._TLanc);
CoolPropDbl rhoV = HEOS.components[0]->ancillaries.rhoV.evaluate(HEOS._TVanc);
CoolPropDbl rhoL = HEOS.components[0].ancillaries.rhoL.evaluate(HEOS._TLanc);
CoolPropDbl rhoV = HEOS.components[0].ancillaries.rhoV.evaluate(HEOS._TVanc);
// Solve for the density
HEOS.SatL->update_TP_guessrho(HEOS._TLanc, HEOS._p, rhoL);
HEOS.SatV->update_TP_guessrho(HEOS._TVanc, HEOS._p, rhoV);
Expand Down Expand Up @@ -709,17 +709,17 @@ void FlashRoutines::PHSU_D_flash(HelmholtzEOSMixtureBackend &HEOS, parameters ot
{
if (HEOS.is_pure_or_pseudopure)
{
CoolPropFluid * component = HEOS.components[0];
CoolPropFluid &component = HEOS.components[0];

shared_ptr<HelmholtzEOSMixtureBackend> Sat;
CoolPropDbl rhoLtriple = component->triple_liquid.rhomolar;
CoolPropDbl rhoVtriple = component->triple_vapor.rhomolar;
CoolPropDbl rhoLtriple = component.triple_liquid.rhomolar;
CoolPropDbl rhoVtriple = component.triple_vapor.rhomolar;
// Check if in the "normal" region
if (HEOS._rhomolar >= rhoVtriple && HEOS._rhomolar <= rhoLtriple)
{
CoolPropDbl yL, yV, value, y_solid;
CoolPropDbl TLtriple = component->triple_liquid.T; ///TODO: separate TL and TV for ppure
CoolPropDbl TVtriple = component->triple_vapor.T;
CoolPropDbl TLtriple = component.triple_liquid.T; ///TODO: separate TL and TV for ppure
CoolPropDbl TVtriple = component.triple_vapor.T;

// First check if solid (below the line connecting the triple point values) - this is an error for now
switch (other)
Expand Down Expand Up @@ -797,10 +797,10 @@ void FlashRoutines::PHSU_D_flash(HelmholtzEOSMixtureBackend &HEOS, parameters ot
}
}
// Check if vapor/solid region below triple point vapor density
else if (HEOS._rhomolar < component->triple_vapor.rhomolar)
else if (HEOS._rhomolar < component.triple_vapor.rhomolar)
{
CoolPropDbl y, value;
CoolPropDbl TVtriple = component->triple_vapor.T; //TODO: separate TL and TV for ppure
CoolPropDbl TVtriple = component.triple_vapor.T; //TODO: separate TL and TV for ppure

// If value is above the value calculated from X(Ttriple, _rhomolar), it is vapor
switch (other)
Expand Down Expand Up @@ -834,7 +834,7 @@ void FlashRoutines::PHSU_D_flash(HelmholtzEOSMixtureBackend &HEOS, parameters ot
else
{
CoolPropDbl y, value;
CoolPropDbl TLtriple = component->pEOS->Ttriple;
CoolPropDbl TLtriple = component.EOS().Ttriple;

// If value is above the value calculated from X(Ttriple, _rhomolar), it is vapor
switch (other)
Expand Down Expand Up @@ -1324,8 +1324,8 @@ void FlashRoutines::HS_flash(HelmholtzEOSMixtureBackend &HEOS)
HEOS.calc_ssat_max();

CoolProp::SimpleState crit = HEOS.get_state("reducing");
CoolProp::SimpleState &tripleL = HEOS.components[0]->triple_liquid,
&tripleV = HEOS.components[0]->triple_vapor;
CoolProp::SimpleState &tripleL = HEOS.components[0].triple_liquid,
&tripleV = HEOS.components[0].triple_vapor;

double first_maxima_in_saturation_entropy;
if (HEOS.ssat_max.exists == SsatSimpleState::SSAT_MAX_DOES_EXIST){
Expand Down
7 changes: 6 additions & 1 deletion src/Backends/Helmholtz/Fluids/FluidLibrary.cpp
Expand Up @@ -23,7 +23,7 @@ JSONFluidLibrary & get_library(void){
return library;
}

CoolPropFluid& get_fluid(const std::string &fluid_string){
CoolPropFluid get_fluid(const std::string &fluid_string){
if (library.is_empty()){ load(); }
return library.get(fluid_string);
}
Expand All @@ -33,4 +33,9 @@ std::string get_fluid_list(void){
return library.get_fluid_list();
};

void set_fluid_enthalpy_entropy_offset(const std::string &fluid, double delta_a1, double delta_a2, const std::string &ref){
if (library.is_empty()){ load(); }
library.set_fluid_enthalpy_entropy_offset(fluid, delta_a1, delta_a2, ref);
}

} /* namespace CoolProp */
32 changes: 23 additions & 9 deletions src/Backends/Helmholtz/Fluids/FluidLibrary.h
Expand Up @@ -380,9 +380,6 @@ class JSONFluidLibrary
{
parse_EOS(*itr,fluid);
}

// Set the EOS pointer to the first EOS
fluid.pEOS = &(fluid.EOSVector[0]);
};

/// Parse the transport properties
Expand Down Expand Up @@ -871,8 +868,8 @@ class JSONFluidLibrary
// Use the method of Chung to approximate the values for epsilon_over_k and sigma_eta
// Chung, T.-H.; Ajlan, M.; Lee, L. L.; Starling, K. E. Generalized Multiparameter Correlation for Nonpolar and Polar Fluid Transport Properties. Ind. Eng. Chem. Res. 1988, 27, 671-679.
// rhoc needs to be in mol/L to yield a sigma in nm,
CoolPropDbl rho_crit_molar = fluid.pEOS->reduce.rhomolar/1000.0;// [mol/m3 to mol/L]
CoolPropDbl Tc = fluid.pEOS->reduce.T;
CoolPropDbl rho_crit_molar = fluid.EOS().reduce.rhomolar/1000.0;// [mol/m3 to mol/L]
CoolPropDbl Tc = fluid.EOS().reduce.T;
fluid.transport.sigma_eta = 0.809/pow(rho_crit_molar, static_cast<CoolPropDbl>(1.0/3.0))/1e9; // 1e9 is to convert from nm to m
fluid.transport.epsilon_over_k = Tc/1.3593; // [K]
}
Expand Down Expand Up @@ -1169,7 +1166,7 @@ class JSONFluidLibrary
/**
@param key Either a CAS number or the name (CAS number should be preferred)
*/
CoolPropFluid& get(const std::string &key)
CoolPropFluid get(const std::string &key)
{
// Try to find it
std::map<std::string, std::size_t>::const_iterator it = string_to_index_map.find(key);
Expand All @@ -1185,7 +1182,7 @@ class JSONFluidLibrary
/**
@param key The index of the fluid in the map
*/
CoolPropFluid& get(std::size_t key)
CoolPropFluid get(std::size_t key)
{
// Try to find it
std::map<std::size_t, CoolPropFluid>::iterator it = fluid_map.find(key);
Expand All @@ -1197,6 +1194,20 @@ class JSONFluidLibrary
throw ValueError(format("key [%d] was not found in JSONFluidLibrary",key));
}
};
void set_fluid_enthalpy_entropy_offset(const std::string &fluid, double delta_a1, double delta_a2, const std::string &ref){
// Try to find it
std::map<std::string, std::size_t>::const_iterator it = string_to_index_map.find(fluid);
if (it != string_to_index_map.end()){
std::map<std::size_t, CoolPropFluid>::iterator it2 = fluid_map.find(it->second);
// If it is found
if (it2 != fluid_map.end()){
it2->second.EOSVector[0].alpha0.EnthalpyEntropyOffset.set(delta_a1, delta_a2, ref);
}
else{
throw ValueError(format("fluid [%s] was not found in JSONFluidLibrary",fluid.c_str()));
}
}
}
/// Return a comma-separated list of fluid names
std::string get_fluid_list(void)
{
Expand All @@ -1210,8 +1221,11 @@ JSONFluidLibrary & get_library(void);
/// Get a comma-separated-list of fluids that are included
std::string get_fluid_list(void);

/// Get the fluid structure returned as a reference
CoolPropFluid& get_fluid(const std::string &fluid_string);
/// Get the fluid structure
CoolPropFluid get_fluid(const std::string &fluid_string);

/// Set the internal enthalpy and entropy offset variables
void set_fluid_enthalpy_entropy_offset(const std::string &fluid, double delta_a1, double delta_a2, const std::string &ref);

} /* namespace CoolProp */
#endif
8 changes: 4 additions & 4 deletions src/Backends/Helmholtz/HelmholtzEOSBackend.h
Expand Up @@ -18,22 +18,22 @@ namespace CoolProp {
class HelmholtzEOSBackend : public HelmholtzEOSMixtureBackend {
public:
HelmholtzEOSBackend();
HelmholtzEOSBackend(CoolPropFluid *pFluid){set_components(std::vector<CoolPropFluid*>(1,pFluid));};
HelmholtzEOSBackend(CoolPropFluid Fluid){set_components(std::vector<CoolPropFluid>(1,Fluid));};
HelmholtzEOSBackend(const std::string &name){
Dictionary dict;
std::vector<double> mole_fractions;
std::vector<CoolPropFluid*> components;
std::vector<CoolPropFluid> components;
if (is_predefined_mixture(name, dict)){
std::vector<std::string> fluids = dict.get_string_vector("fluids");
mole_fractions = dict.get_double_vector("mole_fractions");

components.resize(fluids.size());
for (unsigned int i = 0; i < components.size(); ++i){
components[i] = &(get_library().get(fluids[i]));
components[i] = get_library().get(fluids[i]);
}
}
else{
components.push_back(&(get_library().get(name))); // Until now it's empty
components.push_back(get_library().get(name)); // Until now it's empty
mole_fractions.push_back(1.);
}
// Set the components
Expand Down

0 comments on commit 6b6f2db

Please sign in to comment.