github/CoolProp
1
1
Fork 0
mirror of https://github.com/CoolProp/CoolProp.git synced 2026-02-09 05:15:45 -05:00
Code Issues Packages Projects Releases Wiki Activity

Added code and tests for saturation ancillaries (p, rho', rho'')

Browse Source 
Closes https://github.com/CoolProp/CoolProp/issues/75

Signed-off-by: Ian Bell <ian.h.bell@gmail.com>
This commit is contained in:
Ian Bell
2014-08-23 13:01:56 +02:00
parent 04beedd13a
commit cd8ee9eb52
5 changed files with 138 additions and 9 deletions
Download Patch File Download Diff File Expand all files Collapse all files

3
src/AbstractState.cpp
View File

@@ -330,6 +330,9 @@ double AbstractState::conductivity(void){
double AbstractState::melting_line(int param, int given, double value){
return calc_melting_line(param, given, value);
}
double AbstractState::saturation_ancillary(parameters param, int Q, parameters given, double value){
return calc_saturation_ancillary(param, Q, given, value);
}
double AbstractState::surface_tension(void){
if (!_surface_tension) _surface_tension = calc_surface_tension();
return _surface_tension;

77
src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp
View File

@@ -174,6 +174,64 @@ long double HelmholtzEOSMixtureBackend::calc_molar_mass(void)
}
return summer;
}
long double HelmholtzEOSMixtureBackend::calc_saturation_ancillary(parameters param, int Q, parameters given, double value)
{
if (is_pure_or_pseudopure)
{
if (param == iP && given == iT){
// p = f(T), direct evaluation
switch (Q)
{
case 0:
return components[0]->ancillaries.pL.evaluate(value);
case 1:
return components[0]->ancillaries.pV.evaluate(value);
}
}
else if (param == iT && given == iP){
// T = f(p), inverse evaluation
switch (Q)
{
case 0:
return components[0]->ancillaries.pL.invert(value);
case 1:
return components[0]->ancillaries.pV.invert(value);
}
}
else if (param == iDmolar && given == iT){
// rho = f(T), inverse evaluation
switch (Q)
{
case 0:
return components[0]->ancillaries.rhoL.evaluate(value);
case 1:
return components[0]->ancillaries.rhoV.evaluate(value);
}
}
else if (param == iT && given == iDmolar){
// T = f(rho), inverse evaluation
switch (Q)
{
case 0:
return components[0]->ancillaries.rhoL.invert(value);
case 1:
return components[0]->ancillaries.rhoV.invert(value);
}
}
else{
throw ValueError(format("calc of %s given %s is invalid in calc_saturation_ancillary",
get_parameter_information(param,"short").c_str(),
get_parameter_information(given,"short").c_str()));
}
throw ValueError(format("Q [%d] is invalid in calc_saturation_ancillary", Q));
}
else
{
throw NotImplementedError(format("calc_saturation_ancillary not implemented for mixtures"));
}
}
long double HelmholtzEOSMixtureBackend::calc_melting_line(int param, int given, long double value)
{
if (is_pure_or_pseudopure)
@@ -1624,23 +1682,26 @@ long double HelmholtzEOSMixtureBackend::solver_rho_Tp(long double T, long double
// Calculate a guess value using SRK equation of state
rhomolar_guess = solver_rho_Tp_SRK(T, p, phase);
if (phase == iphase_gas || phase == iphase_supercritical_gas)
// A gas-like phase, ideal gas might not be the perfect model, but probably good enough
if (phase == iphase_gas || phase == iphase_supercritical_gas || iphase_supercritical)
{
if (rhomolar_guess < 0 || !ValidNumber(rhomolar_guess)) // If the guess is bad, probably high temperature, use ideal gas
{
rhomolar_guess = p/(gas_constant()*T);
}
}
else
// It's liquid at subcritical pressure, we can use ancillaries as a backup
else if (phase == iphase_liquid)
{
if (phase == iphase_liquid)
{
long double _rhoLancval = static_cast<long double>(components[0]->ancillaries.rhoL.evaluate(T));
if (!ValidNumber(rhomolar_guess) || rhomolar_guess < _rhoLancval){
rhomolar_guess = _rhoLancval;
}
long double _rhoLancval = static_cast<long double>(components[0]->ancillaries.rhoL.evaluate(T));
if (!ValidNumber(rhomolar_guess) || rhomolar_guess < _rhoLancval){
rhomolar_guess = _rhoLancval;
}
}
else if (phase == iphase_supercritical_liquid){
long double T = static_cast<long double>(components[0]->ancillaries.pL.invert(_p));
rhomolar_guess = components[0]->ancillaries.rhoL.evaluate(T);
}
}
try{

1
src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.h
View File

@@ -53,6 +53,7 @@ public:
bool has_melting_line(){ return is_pure_or_pseudopure && components[0]->ancillaries.melting_line.enabled();};
long double calc_melting_line(int param, int given, long double value);
phases calc_phase(void){return _phase;};
long double calc_saturation_ancillary(parameters param, int Q, parameters given, double value);
const CoolProp::SimpleState &calc_state(const std::string &state);

58
src/Tests/CoolProp-Tests.cpp
View File

@@ -916,6 +916,64 @@ TEST_CASE("Test partial derivatives using PropsSI", "[derivatives]")
}
}
TEST_CASE("Ancillary functions", "[ancillary]")
{
std::vector<std::string> fluids = strsplit(CoolProp::get_global_param_string("fluids_list"),',');
for (std::size_t i = 0; i < fluids.size(); ++i){
shared_ptr<CoolProp::AbstractState> AS(CoolProp::AbstractState::factory("HEOS", fluids[i]));
double Tc = AS->T_critical();
double Tt = AS->Ttriple();
for (double f = 0.1; f < 1; f += 0.2)
{
double T = f*Tc + (1-f)*Tt;
std::ostringstream ss1;
ss1 << "Pressure error < 2% for fluid " << fluids[i] << " at " << T << " K";
SECTION(ss1.str(), "")
{
AS->update(CoolProp::QT_INPUTS, 0, T);
double p_EOS = AS->p();
double p_anc = AS->saturation_ancillary(CoolProp::iP, 0, CoolProp::iT, T);
double err = std::abs(p_EOS-p_anc)/p_anc;
CAPTURE(p_EOS);
CAPTURE(p_anc);
CAPTURE(T);
CHECK(err < 0.02);
}
std::ostringstream ss2;
ss2 << "Liquid density error < 3% for fluid " << fluids[i] << " at " << T << " K";
SECTION(ss2.str(), "")
{
AS->update(CoolProp::QT_INPUTS, 0, T);
double rho_EOS = AS->rhomolar();
double rho_anc = AS->saturation_ancillary(CoolProp::iDmolar, 0, CoolProp::iT, T);
double err = std::abs(rho_EOS-rho_anc)/rho_anc;
CAPTURE(rho_EOS);
CAPTURE(rho_anc);
CAPTURE(T);
CHECK(err < 0.03);
}
std::ostringstream ss3;
ss3 << "Vapor density error < 3% for fluid " << fluids[i] << " at " << T << " K";
SECTION(ss3.str(), "")
{
AS->update(CoolProp::QT_INPUTS, 1, T);
double rho_EOS = AS->rhomolar();
double rho_anc = AS->saturation_ancillary(CoolProp::iDmolar, 1, CoolProp::iT, T);
double err = std::abs(rho_EOS-rho_anc)/rho_anc;
CAPTURE(rho_EOS);
CAPTURE(rho_anc);
CAPTURE(T);
CHECK(err < 0.03);
}
}
}
}
//TEST_CASE("Test that states agree with CoolProp", "[states]")
//{
// std::vector<std::string> fluids = strsplit(CoolProp::get_global_param_string("fluids_list"),',');