Much more progress on PY flash routines. Only a few edge cases that fail, otherwise all is good

Moved from bounded secant to Brent for stability

Signed-off-by: Ian Bell <ian.h.bell@gmail.com>
This commit is contained in:
Ian Bell
2014-08-13 14:26:47 +02:00
parent baf5a88a8e
commit a7ddd80d5a
16 changed files with 118 additions and 113 deletions

View File

@@ -256,17 +256,17 @@
"sat_min_liquid": {
"T": 159.10000000000002,
"T_units": "K",
"p": -1,
"p": 0.0007350277381973985,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 19729.112579315948,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 159.10000000000002,
"T_units": "K",
"p": -1,
"p": 0.0007350470774723961,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 5.556614984609872e-07,
"rhomolar_units": "mol/m^3"
}
},
@@ -488,17 +488,17 @@
"triple_liquid": {
"T": 159.10000000000002,
"T_units": "K",
"p": -1,
"p": 0.0007350277381973985,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 19729.112579315948,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 159.10000000000002,
"T_units": "K",
"p": -1,
"p": 0.0007350470774723961,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 5.556614984609872e-07,
"rhomolar_units": "mol/m^3"
}
},

View File

@@ -257,19 +257,19 @@
"rhomolar_units": "mol/m^3"
},
"sat_min_liquid": {
"T": 53.481100000000005,
"T": 53.48110000000001,
"T_units": "K",
"p": 238.81033157256962,
"p": 238.81033090289185,
"p_units": "Pa",
"rhomolar": 44917.310369072664,
"rhomolar": 44917.31036907279,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 53.481100000000005,
"T": 53.48110000000001,
"T_units": "K",
"p": 238.81033291967333,
"p": 238.81033291967213,
"p_units": "Pa",
"rhomolar": 0.5372485038942366,
"rhomolar": 0.5372485038942338,
"rhomolar_units": "mol/m^3"
}
},
@@ -517,19 +517,19 @@
"rhomolar_units": "mol/m^3"
},
"triple_liquid": {
"T": 53.481100000000005,
"T": 53.48110000000001,
"T_units": "K",
"p": 238.81033157256962,
"p": 238.81033090289185,
"p_units": "Pa",
"rhomolar": 44917.310369072664,
"rhomolar": 44917.31036907279,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 53.481100000000005,
"T": 53.48110000000001,
"T_units": "K",
"p": 238.81033291967333,
"p": 238.81033291967213,
"p_units": "Pa",
"rhomolar": 0.5372485038942366,
"rhomolar": 0.5372485038942338,
"rhomolar_units": "mol/m^3"
}
}

View File

@@ -238,19 +238,19 @@
"rhomolar_units": "mol/m^3"
},
"sat_min_liquid": {
"T": 187.70000000000002,
"T": 187.7,
"T_units": "K",
"p": 23258.85584048954,
"p": 23258.855840772056,
"p_units": "Pa",
"rhomolar": 29116.32684062676,
"rhomolar": 29116.326840626763,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 187.70000000000002,
"T": 187.7,
"T_units": "K",
"p": 23258.85584073585,
"p": 23258.85584073579,
"p_units": "Pa",
"rhomolar": 15.0242285718773,
"rhomolar": 15.024228571877261,
"rhomolar_units": "mol/m^3"
}
},
@@ -373,19 +373,19 @@
"rhomolar_units": "mol/m^3"
},
"triple_liquid": {
"T": 187.70000000000002,
"T": 187.7,
"T_units": "K",
"p": 23258.85584048954,
"p": 23258.855840772056,
"p_units": "Pa",
"rhomolar": 29116.32684062676,
"rhomolar": 29116.326840626763,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 187.70000000000002,
"T": 187.7,
"T_units": "K",
"p": 23258.85584073585,
"p": 23258.85584073579,
"p_units": "Pa",
"rhomolar": 15.0242285718773,
"rhomolar": 15.024228571877261,
"rhomolar_units": "mol/m^3"
}
},

View File

@@ -238,19 +238,19 @@
"rhomolar_units": "mol/m^3"
},
"sat_min_liquid": {
"T": 119.60000000000001,
"T": 119.6,
"T_units": "K",
"p": -1,
"p": 7.626710591531891e-06,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 9370.782450632158,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 119.60000000000001,
"T": 119.6,
"T_units": "K",
"p": -1,
"p": 7.673974446178117e-06,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 7.717106504504615e-09,
"rhomolar_units": "mol/m^3"
}
},
@@ -368,19 +368,19 @@
"rhomolar_units": "mol/m^3"
},
"triple_liquid": {
"T": 119.60000000000001,
"T": 119.6,
"T_units": "K",
"p": -1,
"p": 7.626710591531891e-06,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 9370.782450632158,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 119.60000000000001,
"T": 119.6,
"T_units": "K",
"p": -1,
"p": 7.673974446178117e-06,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 7.717106504504615e-09,
"rhomolar_units": "mol/m^3"
}
}

View File

@@ -256,17 +256,17 @@
"sat_min_liquid": {
"T": 112.65,
"T_units": "K",
"p": -1,
"p": 8.913009802394882e-05,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 10935.8893141883,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 112.65,
"T_units": "K",
"p": -1,
"p": 8.952745179445645e-05,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 9.558513393991252e-08,
"rhomolar_units": "mol/m^3"
}
},
@@ -386,17 +386,17 @@
"triple_liquid": {
"T": 112.65,
"T_units": "K",
"p": -1,
"p": 8.913009802394882e-05,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 10935.8893141883,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 112.65,
"T_units": "K",
"p": -1,
"p": 8.952745179445645e-05,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 9.558513393991252e-08,
"rhomolar_units": "mol/m^3"
}
}

View File

@@ -238,17 +238,17 @@
"rhomolar_units": "mol/m^3"
},
"sat_min_liquid": {
"T": 205.20000000000002,
"T": 205.2,
"T_units": "K",
"p": 0.0004796318407204005,
"p": 0.0004796318407204004,
"p_units": "Pa",
"rhomolar": 3032.1188560391724,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 205.20000000000002,
"T": 205.2,
"T_units": "K",
"p": 0.00047950993795159446,
"p": 0.0004795099379515944,
"p_units": "Pa",
"rhomolar": 2.8105129154408067e-07,
"rhomolar_units": "mol/m^3"
@@ -371,17 +371,17 @@
"rhomolar_units": "mol/m^3"
},
"triple_liquid": {
"T": 205.20000000000002,
"T": 205.2,
"T_units": "K",
"p": 0.0004796318407204005,
"p": 0.0004796318407204004,
"p_units": "Pa",
"rhomolar": 3032.1188560391724,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 205.20000000000002,
"T": 205.2,
"T_units": "K",
"p": 0.00047950993795159446,
"p": 0.0004795099379515944,
"p_units": "Pa",
"rhomolar": 2.8105129154408067e-07,
"rhomolar_units": "mol/m^3"

View File

@@ -238,19 +238,19 @@
"rhomolar_units": "mol/m^3"
},
"sat_min_liquid": {
"T": 187.20000000000002,
"T": 187.2,
"T_units": "K",
"p": 0.0007989708776559431,
"p": 0.000798968839884266,
"p_units": "Pa",
"rhomolar": 3930.823263642575,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 187.20000000000002,
"T": 187.2,
"T_units": "K",
"p": 0.0007991110509185737,
"p": 0.0007991110509185733,
"p_units": "Pa",
"rhomolar": 5.134127178430104e-07,
"rhomolar": 5.134127178430102e-07,
"rhomolar_units": "mol/m^3"
}
},
@@ -371,19 +371,19 @@
"rhomolar_units": "mol/m^3"
},
"triple_liquid": {
"T": 187.20000000000002,
"T": 187.2,
"T_units": "K",
"p": 0.0007989708776559431,
"p": 0.000798968839884266,
"p_units": "Pa",
"rhomolar": 3930.823263642575,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 187.20000000000002,
"T": 187.2,
"T_units": "K",
"p": 0.0007991110509185737,
"p": 0.0007991110509185733,
"p_units": "Pa",
"rhomolar": 5.134127178430104e-07,
"rhomolar": 5.134127178430102e-07,
"rhomolar_units": "mol/m^3"
}
}

View File

@@ -263,17 +263,17 @@
"sat_min_liquid": {
"T": 87.953,
"T_units": "K",
"p": -1,
"p": 0.0007478599372011701,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 18254.562780551252,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 87.953,
"T_units": "K",
"p": -1,
"p": 0.0007471728686973618,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 1.021729100294641e-06,
"rhomolar_units": "mol/m^3"
}
},
@@ -479,17 +479,17 @@
"triple_liquid": {
"T": 87.953,
"T_units": "K",
"p": -1,
"p": 0.0007478599372011701,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 18254.562780551252,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 87.953,
"T_units": "K",
"p": -1,
"p": 0.0007471728686973618,
"p_units": "Pa",
"rhomolar": -1,
"rhomolar": 1.021729100294641e-06,
"rhomolar_units": "mol/m^3"
}
}

View File

@@ -242,19 +242,19 @@
"rhomolar_units": "mol/m^3"
},
"sat_min_liquid": {
"T": 197.70000000000002,
"T": 197.7,
"T_units": "K",
"p": 1660.2201918556232,
"p": 1660.2201919340835,
"p_units": "Pa",
"rhomolar": 25290.03160667635,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 197.70000000000002,
"T": 197.7,
"T_units": "K",
"p": 1660.2201925857255,
"p": 1660.2201925857146,
"p_units": "Pa",
"rhomolar": 1.0119785850078464,
"rhomolar": 1.01197858500784,
"rhomolar_units": "mol/m^3"
}
},
@@ -377,19 +377,19 @@
"rhomolar_units": "mol/m^3"
},
"triple_liquid": {
"T": 197.70000000000002,
"T": 197.7,
"T_units": "K",
"p": 1660.2201918556232,
"p": 1660.2201919340835,
"p_units": "Pa",
"rhomolar": 25290.03160667635,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 197.70000000000002,
"T": 197.7,
"T_units": "K",
"p": 1660.2201925857255,
"p": 1660.2201925857146,
"p_units": "Pa",
"rhomolar": 1.0119785850078464,
"rhomolar": 1.01197858500784,
"rhomolar_units": "mol/m^3"
}
}

View File

@@ -239,17 +239,17 @@
"rhomolar_units": "mol/m^3"
},
"sat_min_liquid": {
"T": 219.70000000000002,
"T": 219.7,
"T_units": "K",
"p": 0.44448529109093327,
"p": 0.44448529109093315,
"p_units": "Pa",
"rhomolar": 6051.9562481855255,
"rhomolar_units": "mol/m^3"
},
"sat_min_vapor": {
"T": 219.70000000000002,
"T": 219.7,
"T_units": "K",
"p": 0.44448541774268197,
"p": 0.4444854177426819,
"p_units": "Pa",
"rhomolar": 0.00024332923985951783,
"rhomolar_units": "mol/m^3"
@@ -369,17 +369,17 @@
"rhomolar_units": "mol/m^3"
},
"triple_liquid": {
"T": 219.70000000000002,
"T": 219.7,
"T_units": "K",
"p": 0.44448529109093327,
"p": 0.44448529109093315,
"p_units": "Pa",
"rhomolar": 6051.9562481855255,
"rhomolar_units": "mol/m^3"
},
"triple_vapor": {
"T": 219.70000000000002,
"T": 219.7,
"T_units": "K",
"p": 0.44448541774268197,
"p": 0.4444854177426819,
"p_units": "Pa",
"rhomolar": 0.00024332923985951783,
"rhomolar_units": "mol/m^3"

View File

@@ -53,12 +53,12 @@
"T_m": 85.39999999999998,
"parts": [
{
"T_0": 85.3,
"T_0": 85.525,
"T_max": 168.63,
"T_min": 85.525,
"a": 718000000.0,
"c": 1.283,
"p_0": 0.0
"p_0": 0.00017200108942021548
}
],
"type": "Simon"

View File

@@ -58,10 +58,11 @@ void FlashRoutines::QT_flash(HelmholtzEOSMixtureBackend &HEOS)
// If you get here, there was no error, all is well
break;
}
catch(std::exception &e){
catch(std::exception &){
if (omega < 1.1*increment){
throw;
}
// else we are going to try again with a smaller omega
}
}
}
@@ -147,7 +148,7 @@ void FlashRoutines::PQ_flash(HelmholtzEOSMixtureBackend &HEOS)
// Check limits
if (!is_in_closed_range(pmin_sat*0.999999, pmax_sat*1.000001, static_cast<long double>(HEOS._p))){
throw ValueError(format("Pressure to PQ_flash [%6g Pa] must be in range [%8g Pa, %8g Pa]",HEOS._p, pmin_sat, pmax_sat));
throw ValueError(format("Pressure to PQ_flash [%6g Pa] must be in range [%8Lg Pa, %8Lg Pa]",HEOS._p, pmin_sat, pmax_sat));
}
// ------------------
// It is a pure fluid
@@ -177,6 +178,7 @@ void FlashRoutines::PQ_flash(HelmholtzEOSMixtureBackend &HEOS)
if (omega < 1.1*increment){
throw;
}
// else we are going to try again with a smaller omega
}
}
}

View File

@@ -101,7 +101,9 @@ double SaturationAncillaryFunction::invert(double value)
std::string errstring;
try{
return Brent(resid,Tmin,Tmax,DBL_EPSILON,1e-12,100,errstring);
// Safe to expand the domain a little bit to lower temperature, absolutely cannot exceed Tmax
// because then you get (negative number)^(double) which is undefined.
return Brent(resid,Tmin-0.01,Tmax,DBL_EPSILON,1e-12,100,errstring);
}
catch(std::exception &e){
return Secant(resid,Tmax, -0.01, 1e-12, 100, errstring);

View File

@@ -1501,6 +1501,7 @@ long double HelmholtzEOSMixtureBackend::solver_rho_Tp(long double T, long double
phase = _phase;
if (rhomolar_guess < 0) // Not provided
{
// 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)
@@ -1514,9 +1515,9 @@ long double HelmholtzEOSMixtureBackend::solver_rho_Tp(long double T, long double
{
if (phase == iphase_liquid)
{
_rhoLanc = components[0]->ancillaries.rhoL.evaluate(T);
if (rhomolar_guess < static_cast<long double>(_rhoLanc)){
rhomolar_guess = static_cast<long double>(_rhoLanc);
long double _rhoLancval = static_cast<long double>(components[0]->ancillaries.rhoL.evaluate(T));
if (!ValidNumber(rhomolar_guess) || rhomolar_guess < _rhoLancval){
rhomolar_guess = _rhoLancval;
}
}
}

View File

@@ -40,14 +40,15 @@ namespace CoolProp {
};
solver_resid resid(HEOS, T, options.rhoL, options.rhoV);
if (!ValidNumber(options.p)){throw ValueError("options.p is not valid in saturation_T_pure_1D_P");};
if (!ValidNumber(options.rhoL)){throw ValueError("options.rhoL is not valid in saturation_T_pure_1D_P");};
if (!ValidNumber(options.rhoV)){throw ValueError("options.rhoV is not valid in saturation_T_pure_1D_P");};
if (!ValidNumber(options.p)){throw ValueError(format("options.p is not valid in saturation_T_pure_1D_P for T = %Lg",T));};
if (!ValidNumber(options.rhoL)){throw ValueError(format("options.rhoL is not valid in saturation_T_pure_1D_P for T = %Lg",T));};
if (!ValidNumber(options.rhoV)){throw ValueError(format("options.rhoV is not valid in saturation_T_pure_1D_P for T = %Lg",T));};
std::string errstr;
long double pmax = std::min(options.p*1.03, static_cast<long double>(HEOS->p_critical()+1e-6));
long double pmin = std::max(options.p*0.97, static_cast<long double>(HEOS->p_triple()-1e-6));
BoundedSecant(resid, options.p, pmin, pmax, sqrt(pmin*pmax), 1e-10, 100, errstr);
Brent(resid, pmin, pmax, LDBL_EPSILON, 1e-8, 100, errstr);
}
void SaturationSolvers::saturation_P_pure_1D_T(HelmholtzEOSMixtureBackend *HEOS, long double p, saturation_PHSU_pure_options &options){
@@ -90,9 +91,9 @@ void SaturationSolvers::saturation_P_pure_1D_T(HelmholtzEOSMixtureBackend *HEOS,
if (!ValidNumber(options.rhoV)){throw ValueError("options.rhoV is not valid in saturation_P_pure_1D_T");};
std::string errstr;
long double Tmax = std::min(options.T + 1, static_cast<long double>(HEOS->T_critical()-1e-6));
long double Tmin = std::max(options.T - 1, static_cast<long double>(HEOS->Ttriple()+1e-6));
BoundedSecant(resid, options.T, Tmin, Tmax, 0.5, 1e-11, 100, errstr);
long double Tmax = std::min(options.T + 2, static_cast<long double>(HEOS->T_critical()-1e-6));
long double Tmin = std::max(options.T - 2, static_cast<long double>(HEOS->Ttriple()+1e-6));
Brent(resid, Tmin, Tmax, LDBL_EPSILON, 1e-11, 100, errstr);
}
void SaturationSolvers::saturation_PHSU_pure(HelmholtzEOSMixtureBackend *HEOS, long double specified_value, saturation_PHSU_pure_options &options)
@@ -478,7 +479,7 @@ void SaturationSolvers::saturation_T_pure(HelmholtzEOSMixtureBackend *HEOS, long
options.rhoL = _options.rhoL;
options.rhoV = _options.rhoV;
options.p = _options.pL;
throw;
SaturationSolvers::saturation_T_pure_1D_P(HEOS, T, options);
}
}
void SaturationSolvers::saturation_T_pure_Akasaka(HelmholtzEOSMixtureBackend *HEOS, long double T, saturation_T_pure_Akasaka_options &options)

View File

@@ -427,7 +427,6 @@ double PropsSI(const std::string &Output, const std::string &Name1, double Prop1
#else
std::cout << "macro is on; error checking disabled in PropsSI" << std::endl;
#endif
// Here is the real code
extract_backend(Ref, backend, fluid);
return _PropsSI(Output, Name1, Prop1, Name2, Prop2, backend, fluid, std::vector<double>());