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

Replaced all fabs() with std::abs()

Browse Source 
Signed-off-by: Ian Bell <ian.h.bell@gmail.com>
This commit is contained in:
Ian Bell
2014-08-23 14:56:07 +02:00
parent b12042bd8c
commit d3261395bf
14 changed files with 133 additions and 133 deletions
Download Patch File Download Diff File Expand all files Collapse all files

8
src/Backends/Helmholtz/FlashRoutines.cpp
View File

@@ -119,7 +119,7 @@ void FlashRoutines::QT_flash(HelmholtzEOSMixtureBackend &HEOS)
rhoLsat = HEOS.solver_rho_Tp(HEOS._T, psatLanc, rhoLanc);
rhoVsat = HEOS.solver_rho_Tp(HEOS._T, psatVanc, rhoVanc);
if (!ValidNumber(rhoLsat) || !ValidNumber(rhoVsat) ||
fabs(rhoLsat/rhoLanc-1) > 0.5 || fabs(rhoVanc/rhoVsat-1) > 0.5)
std::abs(rhoLsat/rhoLanc-1) > 0.5 || std::abs(rhoVanc/rhoVsat-1) > 0.5)
{
throw ValueError("pseudo-pure failed");
}
@@ -524,10 +524,10 @@ void FlashRoutines::HSU_P_flash_singlephase_Newton(HelmholtzEOSMixtureBackend &H
tau -= omega*(B[0][0]*f1+B[0][1]*f2);
delta -= omega*(B[1][0]*f1+B[1][1]*f2);
if (fabs(f1)>fabs(f2))
worst_error=fabs(f1);
if (std::abs(f1) > std::abs(f2))
worst_error = std::abs(f1);
else
worst_error=fabs(f2);
worst_error = std::abs(f2);
if (!ValidNumber(f1) || !ValidNumber(f2))
{

6
src/Backends/Helmholtz/Fluids/FluidLibrary.h
View File

@@ -222,7 +222,7 @@ protected:
long double T0 = cpjson::get_double(contribution, "T0");
// Take the constant term if nonzero and set it as a polyT term
if (fabs(constants[0]) > 1e-14){
if (std::abs(constants[0]) > 1e-14){
std::vector<long double> c(1,constants[0]), t(1,0);
if (EOS.alpha0.CP0PolyT.is_enabled() == true){
EOS.alpha0.CP0PolyT.extend(c,t);
@@ -233,7 +233,7 @@ protected:
}
std::vector<long double> n, c, d, t;
if (fabs(constants[1]) > 1e-14){
if (std::abs(constants[1]) > 1e-14){
// sinh term can be converted by setting a_k = C, b_k = 2*D, c_k = -1, d_k = 1
n.push_back(constants[1]);
t.push_back(-2*constants[2]/Tc);
@@ -241,7 +241,7 @@ protected:
d.push_back(-1);
}
if (fabs(constants[3]) > 1e-14){
if (std::abs(constants[3]) > 1e-14){
// cosh term can be converted by setting a_k = C, b_k = 2*D, c_k = 1, d_k = 1
n.push_back(-constants[3]);
t.push_back(-2*constants[4]/Tc);

10
src/Backends/Helmholtz/TransportRoutines.cpp
View File

@@ -282,14 +282,14 @@ long double TransportRoutines::viscosity_water_hardcoded(HelmholtzEOSMixtureBack
else
{
psi_D=acos(pow(1+powInt(qd*zeta,2),-1.0/2.0));
w=sqrt(fabs((qc*zeta-1)/(qc*zeta+1)))*tan(psi_D/2.0);
w=sqrt(std::abs((qc*zeta-1)/(qc*zeta+1)))*tan(psi_D/2.0);
if (qc*zeta>1){
L=log((1+w)/(1-w));
}
else{
L=2*atan(fabs(w));
L=2*atan(std::abs(w));
}
Y=1.0/12.0*sin(3*psi_D)-1/(4*qc*zeta)*sin(2*psi_D)+1.0/powInt(qc*zeta,2)*(1-5.0/4.0*powInt(qc*zeta,2))*sin(psi_D)-1.0/powInt(qc*zeta,3)*((1-3.0/2.0*powInt(qc*zeta,2))*psi_D-pow(fabs(powInt(qc*zeta,2)-1),3.0/2.0)*L);
Y=1.0/12.0*sin(3*psi_D)-1/(4*qc*zeta)*sin(2*psi_D)+1.0/powInt(qc*zeta,2)*(1-5.0/4.0*powInt(qc*zeta,2))*sin(psi_D)-1.0/powInt(qc*zeta,3)*((1-3.0/2.0*powInt(qc*zeta,2))*psi_D-pow(std::abs(powInt(qc*zeta,2)-1),3.0/2.0)*L);
}
mubar_2=exp(x_mu*Y);
@@ -780,7 +780,7 @@ long double TransportRoutines::conductivity_critical_hardcoded_ammonia(Helmholtz
double pi=3.141592654,a_chi,k_B=1.3806504e-23,X_T,DELTA_lambda,dPdT,eta_B,DELTA_lambda_id,DELTA_lambda_i;
rho = HEOS.keyed_output(CoolProp::iDmass);
t = fabs((T-Tc)/Tc);
t = std::abs((T-Tc)/Tc);
a_chi = a_zeta/0.7;
eta_B = (2.60+1.6*t)*1e-5;
dPdT = (2.18-0.12/exp(17.8*t))*1e5; // [Pa-K]
@@ -831,7 +831,7 @@ long double TransportRoutines::conductivity_hardcoded_helium(HelmholtzEOSMixture
double x0 = 0.392, E1 = 2.8461, E2 = 0.27156, beta = 0.3554, gamma = 1.1743, delta = 4.304, rhoc_crit = 69.158,
Tc = 5.18992, pc = 2.2746e5;
double DeltaT = fabs(1-T/Tc), DeltaRho = fabs(1-rho/rhoc_crit);
double DeltaT = std::abs(1-T/Tc), DeltaRho = std::abs(1-rho/rhoc_crit);
double eta = HEOS.viscosity(); // [Pa-s]
double K_T = HEOS.isothermal_compressibility(), K_Tprime, K_Tbar;
double dpdT = HEOS.first_partial_deriv(CoolProp::iP, CoolProp::iT, CoolProp::iDmolar);

12
src/Backends/Helmholtz/VLERoutines.cpp
View File

@@ -524,7 +524,7 @@ void SaturationSolvers::saturation_D_pure(HelmholtzEOSMixtureBackend *HEOS, long
}
while (error > 1e-9);
long double p_error_limit = 1e-3;
if (fabs(p_error) > p_error_limit){
if (std::abs(p_error) > p_error_limit){
throw SolutionError(format("saturation_D_pure solver abs error on p [%Lg] > limit [%Lg]", p_error, p_error_limit));
}
}
@@ -673,16 +673,16 @@ void SaturationSolvers::saturation_T_pure_Akasaka(HelmholtzEOSMixtureBackend *HE
throw SolutionError(format("Akasaka solver did not converge after 100 iterations"));
}
}
while (error > 1e-10 && fabs(stepL) > 10*DBL_EPSILON*fabs(stepL) && fabs(stepV) > 10*DBL_EPSILON*fabs(stepV));
while (error > 1e-10 && std::abs(stepL) > 10*DBL_EPSILON*std::abs(stepL) && std::abs(stepV) > 10*DBL_EPSILON*std::abs(stepV));
long double p_error_limit = 1e-3;
long double p_error = (PL - PV)/PL;
if (fabs(p_error) > p_error_limit){
if (std::abs(p_error) > p_error_limit){
options.pL = PL;
options.pV = PV;
options.rhoL = rhoL;
options.rhoV = rhoV;
throw SolutionError(format("saturation_T_pure_Akasaka solver abs error on p [%g] > limit [%g]", fabs(p_error), p_error_limit));
throw SolutionError(format("saturation_T_pure_Akasaka solver abs error on p [%g] > limit [%g]", std::abs(p_error), p_error_limit));
}
}
@@ -774,7 +774,7 @@ void SaturationSolvers::successive_substitution(HelmholtzEOSMixtureBackend &HEOS
throw ValueError(format("saturation_p was unable to reach a solution within 50 iterations"));
}
}
while(fabs(f) > 1e-12 || iter < options.Nstep_max);
while(std::abs(f) > 1e-12 || iter < options.Nstep_max);
HEOS.SatL->update_TP_guessrho(T, p, rhomolar_liq);
HEOS.SatV->update_TP_guessrho(T, p, rhomolar_vap);
@@ -910,7 +910,7 @@ void SaturationSolvers::newton_raphson_VLE_GV::call(HelmholtzEOSMixtureBackend *
IO.T = exp(log(IO.T) + v[N]);
IO.rhomolar_liq = exp(log(IO.rhomolar_liq) + v[N+1]);
if (fabs(IO.T) > 1e6)
if (std::abs(IO.T) > 1e6)
{
/*std::cout << "J = " << vec_to_string(J,"%16.15g");
std::cout << "nr = " << vec_to_string(r,"%16.15g");*/

8
src/Backends/REFPROP/REFPROPMixtureBackend.cpp
View File

@@ -1305,7 +1305,7 @@ TEST_CASE("Check REFPROP and CoolProp values agree","[REFPROP]")
CAPTURE(rho_RP);
double DH = (rho_RP-rho_CP)/rho_RP;
CHECK(fabs(DH) < 0.005);
CHECK(std::abs(DH) < 0.005);
}
}
SECTION("Saturation specific heats agree within 0.5% at T/Tc = 0.9")
@@ -1336,7 +1336,7 @@ TEST_CASE("Check REFPROP and CoolProp values agree","[REFPROP]")
CAPTURE(0.9*Tr);
double Dcp = (cp_RP-cp_CP)/cp_RP;
CHECK(fabs(Dcp) < 0.005);
CHECK(std::abs(Dcp) < 0.005);
}
}
SECTION("Enthalpy and entropy reference state")
@@ -1371,8 +1371,8 @@ TEST_CASE("Check REFPROP and CoolProp values agree","[REFPROP]")
double DH = (S1->hmass()-S2->hmass());
double DS = (S1->smass()-S2->smass());
CHECK(fabs(DH/h_RP) < 0.001);
CHECK(fabs(DS/s_RP) < 0.001);
CHECK(std::abs(DH/h_RP) < 0.001);
CHECK(std::abs(DS/s_RP) < 0.001);
}
}
}

2
src/CoolPropTools.cpp
View File

@@ -175,7 +175,7 @@ void solve_cubic(double a, double b, double c, double d, int &N, double &x0, dou
double t0;
if (4*p*p*p+27*q*q>0 && p<0)
{
t0 = -2.0*fabs(q)/q*sqrt(-p/3.0)*cosh(1.0/3.0*acosh(-3.0*fabs(q)/(2.0*p)*sqrt(-3.0/p)));
t0 = -2.0*std::abs(q)/q*sqrt(-p/3.0)*cosh(1.0/3.0*acosh(-3.0*std::abs(q)/(2.0*p)*sqrt(-3.0/p)));
}
else
{

38
src/Helmholtz.cpp
View File

@@ -305,7 +305,7 @@ void ResidualHelmholtzNonAnalytic::all(const long double &tau, const long double
long double dDELTA3_dDelta2_dTau = 2.0*Ai*(betai-1)/(betai*betai)*pow(pow(delta-1,2),1/(2*betai)-1.0);
long double dDELTAbi_dDelta, dDELTA2_dDelta2, dDELTAbi2_dDelta2, dDELTAbi3_dDelta3, dDELTA3_dDelta3;
if (fabs(delta-1) < 10*DBL_EPSILON){
if (std::abs(delta-1) < 10*DBL_EPSILON){
dDELTAbi_dDelta = 0;
dDELTA2_dDelta2 = 0;
dDELTA3_dDelta3 = 0;
@@ -382,7 +382,7 @@ long double ResidualHelmholtzNonAnalytic::dDelta(const long double &tau, const l
long double dDELTA_dDelta=(delta-1.0)*(Ai*theta*2.0/betai*pow(pow(delta-1.0,2),1.0/(2.0*betai)-1.0)+2.0*Bi*ai*pow(pow(delta-1.0,2),ai-1.0));
// At critical point, DELTA is 0, and 1/0^n is undefined
if (fabs(DELTA) < 10*DBL_EPSILON)
if (std::abs(DELTA) < 10*DBL_EPSILON)
{
dDELTAbi_dDelta = 0;
}
@@ -432,7 +432,7 @@ long double ResidualHelmholtzNonAnalytic::dDelta2(const long double &tau, const
long double dPSI2_dDelta2=(2.0*Ci*pow(delta-1.0,2)-1.0)*2.0*Ci*PSI;
if (fabs(delta-1) < 10*DBL_EPSILON){
if (std::abs(delta-1) < 10*DBL_EPSILON){
dDELTA2_dDelta2 = 0;
dDELTAbi2_dDelta2 = 0;
}
@@ -442,7 +442,7 @@ long double ResidualHelmholtzNonAnalytic::dDelta2(const long double &tau, const
}
// At critical point, DELTA is 0, and 1/0^n is undefined
if (fabs(DELTA) < 10*DBL_EPSILON)
if (std::abs(DELTA) < 10*DBL_EPSILON)
{
dDELTAbi_dDelta = 0;
}
@@ -476,7 +476,7 @@ long double ResidualHelmholtzNonAnalytic::dDelta_dTau(const long double &tau, co
long double dDELTAbi_dTau=-2.0*theta*bi*pow(DELTA,bi-1.0);
// At critical point, DELTA is 0, and 1/0^n is undefined
if (fabs(DELTA) < 10*DBL_EPSILON)
if (std::abs(DELTA) < 10*DBL_EPSILON)
{
dDELTAbi_dDelta = 0;
}
@@ -537,7 +537,7 @@ long double ResidualHelmholtzNonAnalytic::dDelta3(const long double &tau, const
long double dDELTAbi3_dDelta3 = bi*(pow(DELTA,bi-1)*dDELTA3_dDelta3+dDELTA2_dDelta2*(bi-1)*pow(DELTA,bi-2)*dDELTA_dDelta+(bi-1)*(pow(DELTA,bi-2)*2*dDELTA_dDelta*dDELTA2_dDelta2+pow(dDELTA_dDelta,2)*(bi-2)*pow(DELTA,bi-3)*dDELTA_dDelta));
// At critical point, DELTA is 0, and 1/0^n is undefined
if (fabs(DELTA) < 10*DBL_EPSILON)
if (std::abs(DELTA) < 10*DBL_EPSILON)
{
dDELTAbi_dDelta = 0;
}
@@ -580,7 +580,7 @@ long double ResidualHelmholtzNonAnalytic::dDelta_dTau2(const long double &tau, c
long double dDELTAbi3_dDelta_dTau2 = 2*bi*(bi-1)*pow(DELTA,bi-2)*dDELTA_dDelta+4*pow(theta,2)*bi*(bi-1)*(bi-2)*pow(DELTA,bi-3)*dDELTA_dDelta+8*theta*bi*(bi-1)*pow(DELTA,bi-2)*dtheta_dDelta;
// At critical point, DELTA is 0, and 1/0^n is undefined
if (fabs(DELTA) < 10*DBL_EPSILON)
if (std::abs(DELTA) < 10*DBL_EPSILON)
{
dDELTAbi_dDelta = 0;
}
@@ -620,7 +620,7 @@ long double ResidualHelmholtzNonAnalytic::dDelta2_dTau(const long double &tau, c
long double dDELTAbi2_dDelta_dTau=-Ai*bi*2.0/betai*pow(DELTA,bi-1.0)*(delta-1.0)*pow(pow(delta-1.0,2),1.0/(2.0*betai)-1.0)-2.0*theta*bi*(bi-1.0)*pow(DELTA,bi-2.0)*dDELTA_dDelta;
// At critical point, DELTA is 0, and 1/0^n is undefined
if (fabs(DELTA) < 10*DBL_EPSILON)
if (std::abs(DELTA) < 10*DBL_EPSILON)
{
dDELTAbi_dDelta = 0;
}
@@ -930,11 +930,11 @@ long double IdealHelmholtzCP0PolyT::base(const long double &tau, const long doub
{
double sum=0;
for (std::size_t i = 0; i < N; ++i){
if (fabs(t[i])<10*DBL_EPSILON)
if (std::abs(t[i])<10*DBL_EPSILON)
{
sum += c[i]-c[i]*tau/tau0+c[i]*log(tau/tau0);
}
else if (fabs(t[i]+1) < 10*DBL_EPSILON)
else if (std::abs(t[i]+1) < 10*DBL_EPSILON)
{
sum += c[i]*tau/Tc*log(tau0/tau)+c[i]/Tc*(tau-tau0);
}
@@ -949,11 +949,11 @@ long double IdealHelmholtzCP0PolyT::dTau(const long double &tau, const long doub
{
double sum=0;
for (std::size_t i = 0; i < N; ++i){
if (fabs(t[i])<10*DBL_EPSILON)
if (std::abs(t[i])<10*DBL_EPSILON)
{
sum += c[i]/tau-c[i]/tau0;
}
else if (fabs(t[i]+1) < 10*DBL_EPSILON)
else if (std::abs(t[i]+1) < 10*DBL_EPSILON)
{
sum += c[i]/Tc*log(tau0/tau);
}
@@ -968,11 +968,11 @@ long double IdealHelmholtzCP0PolyT::dTau2(const long double &tau, const long dou
{
double sum=0;
for (std::size_t i = 0; i < N; ++i){
if (fabs(t[i])<10*DBL_EPSILON)
if (std::abs(t[i])<10*DBL_EPSILON)
{
sum += -c[i]/(tau*tau);
}
else if (fabs(t[i]+1) < 10*DBL_EPSILON)
else if (std::abs(t[i]+1) < 10*DBL_EPSILON)
{
sum += -c[i]/(tau*Tc);
}
@@ -987,11 +987,11 @@ long double IdealHelmholtzCP0PolyT::dTau3(const long double &tau, const long dou
{
double sum=0;
for (std::size_t i = 0; i < N; ++i){
if (fabs(t[i])<10*DBL_EPSILON)
if (std::abs(t[i])<10*DBL_EPSILON)
{
sum += 2*c[i]/(tau*tau*tau);
}
else if (fabs(t[i]+1) < 10*DBL_EPSILON)
else if (std::abs(t[i]+1) < 10*DBL_EPSILON)
{
sum += c[i]/(tau*tau*Tc);
}
@@ -1272,11 +1272,11 @@ public:
};
double err(double v1, double v2)
{
if (fabs(v2) > 1e-15){
return fabs((v1-v2)/v2);
if (std::abs(v2) > 1e-15){
return std::abs((v1-v2)/v2);
}
else{
return fabs(v1-v2);
return std::abs(v1-v2);
}
}
};

96
src/HumidAirProp.cpp
View File

@@ -184,7 +184,7 @@ static double Secant_HAProps_T(const std::string &OutputName, const std::string
int iter=1;
std::string sT = "T";
while ((iter<=3 || (fabs(f)>eps && fabs(change)>1e-10)) && iter<100)
while ((iter<=3 || (std::abs(f)>eps && std::abs(change)>1e-10)) && iter<100)
{
if (iter==1){x1=T_guess; T=x1;}
if (iter==2){x2=T_guess+0.001; T=x2;}
@@ -209,7 +209,7 @@ static double Secant_HAProps_W(const std::string &OutputName, const std::string
double x1=0,x2=0,x3=0,y1=0,y2=0,eps=1e-8,f=999,W=0.0001;
int iter=1;
while ((iter<=3 || fabs(f)>eps) && iter<100)
while ((iter<=3 || std::abs(f)>eps) && iter<100)
{
if (iter == 1){x1 = W_guess; W = x1;}
if (iter == 2){x2 = W_guess+0.001; W = x2;}
@@ -522,7 +522,7 @@ double f_factor(double T, double p)
{
y2=LHS-RHS;
x3=x2-y2/(y2-y1)*(x2-x1);
change=fabs(y2/(y2-y1)*(x2-x1));
change=std::abs(y2/(y2-y1)*(x2-x1));
y1=y2; x1=x2; x2=x3;
}
iter=iter+1;
@@ -619,7 +619,7 @@ double MolarVolume(double T, double p, double psi_w)
Cm=C_m(T,psi_w);
iter=1; eps=1e-11; resid=999;
while ((iter<=3 || fabs(resid)>eps) && iter<100)
while ((iter<=3 || std::abs(resid)>eps) && iter<100)
{
if (iter==1){x1=v_bar0; v_bar=x1;}
if (iter==2){x2=v_bar0+0.000001; v_bar=x2;}
@@ -788,7 +788,7 @@ double MolarEntropy(double T, double p, double psi_w, double v_bar)
vbar_a_guess = R_bar_Lem*T/p; //[m^3/mol] since p in [Pa]
while ((iter<=3 || fabs(f)>eps) && iter<100)
while ((iter<=3 || std::abs(f)>eps) && iter<100)
{
if (iter==1){x1=vbar_a_guess; vbar_a=x1;}
if (iter==2){x2=vbar_a_guess+0.001; vbar_a=x2;}
@@ -853,7 +853,7 @@ double DewpointTemperature(double T, double p, double psi_w)
// p_w_s = p_w, and get guess for T from saturation temperature
iter=1; eps=1e-8; resid=999;
while ((iter<=3 || fabs(resid)>eps) && iter<100)
while ((iter<=3 || std::abs(resid)>eps) && iter<100)
{
if (iter==1){x1 = T0; Tdp=x1;}
if (iter==2){x2 = x1 + 0.1; Tdp=x2;}
@@ -1303,7 +1303,7 @@ double HAPropsSI(const std::string &OutputName, const std::string &Input1Name, d
try{
T = Secant_HAProps_T(SecondaryInputName,(char *)"P",p,MainInputName,MainInputValue,SecondaryInputValue,T_guess);
double val = HAPropsSI(SecondaryInputName,(char *)"T",T,(char *)"P",p,MainInputName,MainInputValue);
if (!ValidNumber(T) || !ValidNumber(val) || !(T_min < T && T < T_max) || fabs(val-SecondaryInputValue)>1e-6)
if (!ValidNumber(T) || !ValidNumber(val) || !(T_min < T && T < T_max) || std::abs(val-SecondaryInputValue)>1e-6)
{
throw CoolProp::ValueError();
}
@@ -1661,73 +1661,73 @@ TEST_CASE("Check HA Virials from Table A.2.1","[RP1485]")
{
SECTION("B_aa")
{
CHECK(fabs(HumidAir::B_Air(-60+273.15)/(-33.065/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::B_Air(0+273.15)/(-13.562/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::B_Air(200+273.15)/(11.905/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::B_Air(350+273.15)/(18.949/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::B_Air(-60+273.15)/(-33.065/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::B_Air(0+273.15)/(-13.562/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::B_Air(200+273.15)/(11.905/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::B_Air(350+273.15)/(18.949/1e6)-1) < 1e-3);
}
SECTION("B_ww")
{
CHECK(fabs(HumidAir::B_Water(-60+273.15)/(-11174/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::B_Water(0+273.15)/(-2025.6/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::B_Water(200+273.15)/(-200.52/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::B_Water(350+273.15)/(-89.888/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::B_Water(-60+273.15)/(-11174/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::B_Water(0+273.15)/(-2025.6/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::B_Water(200+273.15)/(-200.52/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::B_Water(350+273.15)/(-89.888/1e6)-1) < 1e-3);
}
SECTION("B_aw")
{
CHECK(fabs(HumidAir::_B_aw(-60+273.15)/(-68.306/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::_B_aw(0+273.15)/(-38.074/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::_B_aw(200+273.15)/(-2.0472/1e6)-1) < 1e-3);
CHECK(fabs(HumidAir::_B_aw(350+273.15)/(7.5200/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::_B_aw(-60+273.15)/(-68.306/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::_B_aw(0+273.15)/(-38.074/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::_B_aw(200+273.15)/(-2.0472/1e6)-1) < 1e-3);
CHECK(std::abs(HumidAir::_B_aw(350+273.15)/(7.5200/1e6)-1) < 1e-3);
}
SECTION("C_aaa")
{
CHECK(fabs(HumidAir::C_Air(-60+273.15)/(2177.9/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::C_Air(0+273.15)/(1893.1/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::C_Air(200+273.15)/(1551.2/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::C_Air(350+273.15)/(1464.7/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::C_Air(-60+273.15)/(2177.9/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::C_Air(0+273.15)/(1893.1/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::C_Air(200+273.15)/(1551.2/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::C_Air(350+273.15)/(1464.7/1e12)-1) < 1e-3);
}
SECTION("C_www")
{
CHECK(fabs(HumidAir::C_Water(-60+273.15)/(-1.5162999202e-04)-1) < 1e-3); // Relaxed criterion for this parameter
CHECK(fabs(HumidAir::C_Water(0+273.15)/(-10981960/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::C_Water(200+273.15)/(-0.00000003713759442)-1) < 1e-3);
CHECK(fabs(HumidAir::C_Water(350+273.15)/(-0.000000001198914198)-1) < 1e-3);
CHECK(std::abs(HumidAir::C_Water(-60+273.15)/(-1.5162999202e-04)-1) < 1e-3); // Relaxed criterion for this parameter
CHECK(std::abs(HumidAir::C_Water(0+273.15)/(-10981960/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::C_Water(200+273.15)/(-0.00000003713759442)-1) < 1e-3);
CHECK(std::abs(HumidAir::C_Water(350+273.15)/(-0.000000001198914198)-1) < 1e-3);
}
SECTION("C_aaw")
{
CHECK(fabs(HumidAir::_C_aaw(-60+273.15)/(1027.3/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::_C_aaw(0+273.15)/(861.02/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::_C_aaw(200+273.15)/(627.15/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::_C_aaw(350+273.15)/(583.79/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::_C_aaw(-60+273.15)/(1027.3/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::_C_aaw(0+273.15)/(861.02/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::_C_aaw(200+273.15)/(627.15/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::_C_aaw(350+273.15)/(583.79/1e12)-1) < 1e-3);
}
SECTION("C_aww")
{
CHECK(fabs(HumidAir::_C_aww(-60+273.15)/(-1821432/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::_C_aww(0+273.15)/(-224234/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::_C_aww(200+273.15)/(-8436.5/1e12)-1) < 1e-3);
CHECK(fabs(HumidAir::_C_aww(350+273.15)/(-2486.9/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::_C_aww(-60+273.15)/(-1821432/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::_C_aww(0+273.15)/(-224234/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::_C_aww(200+273.15)/(-8436.5/1e12)-1) < 1e-3);
CHECK(std::abs(HumidAir::_C_aww(350+273.15)/(-2486.9/1e12)-1) < 1e-3);
}
}
TEST_CASE("Enhancement factor from Table A.3","[RP1485]")
{
CHECK(fabs(HumidAir::f_factor(-60+273.15,101325)/(1.00708)-1) < 1e-3);
CHECK(fabs(HumidAir::f_factor( 80+273.15,101325)/(1.00573)-1) < 1e-3);
CHECK(fabs(HumidAir::f_factor(-60+273.15,10000e3)/(2.23918)-1) < 1e-3);
CHECK(fabs(HumidAir::f_factor(300+273.15,10000e3)/(1.04804)-1) < 1e-3);
CHECK(std::abs(HumidAir::f_factor(-60+273.15,101325)/(1.00708)-1) < 1e-3);
CHECK(std::abs(HumidAir::f_factor( 80+273.15,101325)/(1.00573)-1) < 1e-3);
CHECK(std::abs(HumidAir::f_factor(-60+273.15,10000e3)/(2.23918)-1) < 1e-3);
CHECK(std::abs(HumidAir::f_factor(300+273.15,10000e3)/(1.04804)-1) < 1e-3);
}
TEST_CASE("Isothermal compressibility from Table A.5","[RP1485]")
{
CHECK(fabs(HumidAir::isothermal_compressibility(-60+273.15,101325)/(0.10771e-9)-1) < 1e-3);
CHECK(fabs(HumidAir::isothermal_compressibility( 80+273.15,101325)/(0.46009e-9)-1) < 1e-2); // Relaxed criterion for this parameter
CHECK(fabs(HumidAir::isothermal_compressibility(-60+273.15,10000e3)/(0.10701e-9)-1) < 1e-3);
CHECK(fabs(HumidAir::isothermal_compressibility(300+273.15,10000e3)/(3.05896e-9)-1) < 1e-3);
CHECK(std::abs(HumidAir::isothermal_compressibility(-60+273.15,101325)/(0.10771e-9)-1) < 1e-3);
CHECK(std::abs(HumidAir::isothermal_compressibility( 80+273.15,101325)/(0.46009e-9)-1) < 1e-2); // Relaxed criterion for this parameter
CHECK(std::abs(HumidAir::isothermal_compressibility(-60+273.15,10000e3)/(0.10701e-9)-1) < 1e-3);
CHECK(std::abs(HumidAir::isothermal_compressibility(300+273.15,10000e3)/(3.05896e-9)-1) < 1e-3);
}
TEST_CASE("Henry constant from Table A.6","[RP1485]")
{
CHECK(fabs(HumidAir::HenryConstant(0.010001+273.15)/(0.22600e-9)-1) < 1e-3);
CHECK(fabs(HumidAir::HenryConstant(300+273.15)/(0.58389e-9)-1) < 1e-3);
CHECK(std::abs(HumidAir::HenryConstant(0.010001+273.15)/(0.22600e-9)-1) < 1e-3);
CHECK(std::abs(HumidAir::HenryConstant(300+273.15)/(0.58389e-9)-1) < 1e-3);
}
// A structure to hold the values for one call to HAProps
@@ -1815,7 +1815,7 @@ TEST_CASE_METHOD(HAPropsConsistencyFixture, "ASHRAE RP1485 Tables", "[RP1485]")
CAPTURE(expected);
std::string errmsg = CoolProp::get_global_param_string("errstring");
CAPTURE(errmsg);
CHECK(fabs(actual/expected-1) < 0.01);
CHECK(std::abs(actual/expected-1) < 0.01);
}
}
SECTION("Table A.11")
@@ -1835,7 +1835,7 @@ TEST_CASE_METHOD(HAPropsConsistencyFixture, "ASHRAE RP1485 Tables", "[RP1485]")
CAPTURE(expected);
std::string errmsg = CoolProp::get_global_param_string("errstring");
CAPTURE(errmsg);
CHECK(fabs(actual/expected-1) < 0.01);
CHECK(std::abs(actual/expected-1) < 0.01);
}
}
SECTION("Table A.12")
@@ -1855,7 +1855,7 @@ TEST_CASE_METHOD(HAPropsConsistencyFixture, "ASHRAE RP1485 Tables", "[RP1485]")
CAPTURE(expected);
std::string errmsg = CoolProp::get_global_param_string("errstring");
CAPTURE(errmsg);
CHECK(fabs(actual/expected-1) < 0.01);
CHECK(std::abs(actual/expected-1) < 0.01);
}
}

8
src/MatrixMath.cpp
View File

@@ -84,26 +84,26 @@ TEST_CASE("Internal consistency checks and example use cases for MatrixMath.h","
Eigen::MatrixXd matrix = CoolProp::vec_to_eigen(vec2D);
for (size_t i = 0; i < matrix.cols(); ++i) {
for (size_t j = 0; j < matrix.rows(); ++j) {
CHECK( fabs(matrix(j,i)-vec2D[j][i]) <= 1e-10 );
CHECK( std::abs(matrix(j,i)-vec2D[j][i]) <= 1e-10 );
}
}
vec2D = CoolProp::eigen_to_vec(matrix);
for (size_t i = 0; i < matrix.cols(); ++i) {
for (size_t j = 0; j < matrix.rows(); ++j) {
CHECK( fabs(matrix(j,i)-vec2D[j][i]) <= 1e-10 );
CHECK( std::abs(matrix(j,i)-vec2D[j][i]) <= 1e-10 );
}
}
std::vector<double> vec1D(cHeat);
matrix = CoolProp::vec_to_eigen(vec1D);
for (size_t i = 0; i < matrix.cols(); ++i) {
for (size_t j = 0; j < matrix.rows(); ++j) {
CHECK( fabs(matrix(j,i)-vec1D[j]) <= 1e-10 );
CHECK( std::abs(matrix(j,i)-vec1D[j]) <= 1e-10 );
}
}
vec1D = CoolProp::eigen_to_vec1D(matrix);
for (size_t i = 0; i < matrix.cols(); ++i) {
for (size_t j = 0; j < matrix.rows(); ++j) {
CHECK( fabs(matrix(j,i)-vec1D[j]) <= 1e-10 );
CHECK( std::abs(matrix(j,i)-vec1D[j]) <= 1e-10 );
}
}
}

10
src/PolyMath.cpp
View File

@@ -457,7 +457,7 @@ double Polynomial2DFrac::evaluate(const Eigen::MatrixXd &coefficients, const dou
if ( (r!=1) && (c!=1) ) {
throw ValueError(format("%s (%d): You have a 2D coefficient matrix (%d,%d), please use the 2D functions. ",__FILE__,__LINE__,coefficients.rows(),coefficients.cols()));
}
if ( (firstExponent<0) && (fabs(x_in-x_base)<DBL_EPSILON)) {
if ( (firstExponent<0) && (std::abs(x_in-x_base)<DBL_EPSILON)) {
throw ValueError(format("%s (%d): A fraction cannot be evaluated with zero as denominator, x_in-x_base=%f ",__FILE__,__LINE__,x_in-x_base));
}
@@ -500,10 +500,10 @@ double Polynomial2DFrac::evaluate(const Eigen::MatrixXd &coefficients, const dou
/// @param x_base double value that represents the base value for a centred fit in the 1st dimension
/// @param y_base double value that represents the base value for a centred fit in the 2nd dimension
double Polynomial2DFrac::evaluate(const Eigen::MatrixXd &coefficients, const double &x_in, const double &y_in, const int &x_exp, const int &y_exp, const double &x_base, const double &y_base){
if ( (x_exp<0) && (fabs(x_in-x_base)<DBL_EPSILON)) {
if ( (x_exp<0) && (std::abs(x_in-x_base)<DBL_EPSILON)) {
throw ValueError(format("%s (%d): A fraction cannot be evaluated with zero as denominator, x_in-x_base=%f ",__FILE__,__LINE__,x_in-x_base));
}
if ( (y_exp<0) && (fabs(y_in-y_base)<DBL_EPSILON)) {
if ( (y_exp<0) && (std::abs(y_in-y_base)<DBL_EPSILON)) {
throw ValueError(format("%s (%d): A fraction cannot be evaluated with zero as denominator, y_in-y_base=%f ",__FILE__,__LINE__,y_in-y_base));
}
@@ -633,7 +633,7 @@ double Polynomial2DFrac::integral(const Eigen::MatrixXd &coefficients, const dou
size_t c = newCoefficients.cols();
if (int_exp==-1) {
if (fabs(int_base)<DBL_EPSILON){
if (std::abs(int_base)<DBL_EPSILON){
Eigen::MatrixXd tmpCoefficients = newCoefficients.row(0) * log(int_val-int_base);
newCoefficients = integrateCoeffs(newCoefficients.block(1,0,r-1,c), 0, 1);
newCoefficients.row(0) = tmpCoefficients;
@@ -1044,7 +1044,7 @@ TEST_CASE("Internal consistency checks and example use cases for PolyMath.cpp","
double x = 0.3, y = 255.3, val1, val2, val3, val4;
//CHECK( fabs( polyInt.derivative(x,y,0)-poly2D.evaluate(x,y) ) <= 1e-10 );
//CHECK( std::abs( polyInt.derivative(x,y,0)-poly2D.evaluate(x,y) ) <= 1e-10 );
std::string tmpStr;

40
src/Solvers.cpp
View File

@@ -30,7 +30,7 @@ std::vector<double> NDNewtonRaphson_Jacobian(FuncWrapperND *f, std::vector<doubl
std::vector<double> f0,v,negative_f0;
std::vector<std::vector<double> > J;
double error = 999;
while (iter==0 || fabs(error)>tol){
while (iter==0 || std::abs(error)>tol){
f0 = f->call(x0);
J = f->Jacobian(x0);
@@ -68,7 +68,7 @@ double Newton(FuncWrapper1D* f, double x0, double ftol, int maxiter, std::string
int iter=1;
errstring.clear();
x = x0;
while (iter < 2 || fabs(fval) > ftol)
while (iter < 2 || std::abs(fval) > ftol)
{
fval = f->call(x);
dx = -fval/f->deriv(x);
@@ -79,7 +79,7 @@ double Newton(FuncWrapper1D* f, double x0, double ftol, int maxiter, std::string
x += dx;
if (fabs(dx/x) < 10*DBL_EPSILON)
if (std::abs(dx/x) < 10*DBL_EPSILON)
{
return x;
}
@@ -116,8 +116,8 @@ double Secant(FuncWrapper1D* f, double x0, double dx, double tol, int maxiter, s
int iter=1;
errstring = "";
if (fabs(dx)==0){ errstring="dx cannot be zero"; return _HUGE;}
while (iter<=2 || fabs(fval)>tol)
if (std::abs(dx)==0){ errstring="dx cannot be zero"; return _HUGE;}
while (iter<=2 || std::abs(fval)>tol)
{
if (iter==1){x1=x0; x=x1;}
if (iter==2){x2=x0+dx; x=x2;}
@@ -137,9 +137,9 @@ double Secant(FuncWrapper1D* f, double x0, double dx, double tol, int maxiter, s
if (iter>1)
{
double deltax = x2-x1;
if (fabs(deltax)<1e-14)
if (std::abs(deltax)<1e-14)
{
if (fabs(fval) < tol*10)
if (std::abs(fval) < tol*10)
{
return x;
}
@@ -182,8 +182,8 @@ double BoundedSecant(FuncWrapper1D* f, double x0, double xmin, double xmax, doub
int iter=1;
errstring = "";
if (fabs(dx)==0){ errstring = "dx cannot be zero"; return _HUGE;}
while (iter<=3 || fabs(fval)>tol)
if (std::abs(dx)==0){ errstring = "dx cannot be zero"; return _HUGE;}
while (iter<=3 || std::abs(fval)>tol)
{
if (iter==1){x1=x0; x=x1;}
else if (iter==2){x2=x0+dx; x=x2;}
@@ -241,11 +241,11 @@ double Brent(FuncWrapper1D* f, double a, double b, double macheps, double t, int
fb = f->call(b);
// If one of the boundaries is to within tolerance, just stop
if (fabs(fb) < t) { return b;}
if (std::abs(fb) < t) { return b;}
if (!ValidNumber(fb)){
throw ValueError(format("Brent's method f(b) is NAN for b = %g",b).c_str());
}
if (fabs(fa) < t) { return a;}
if (std::abs(fa) < t) { return a;}
if (!ValidNumber(fa)){
throw ValueError(format("Brent's method f(a) is NAN for a = %g",a).c_str());
}
@@ -256,7 +256,7 @@ double Brent(FuncWrapper1D* f, double a, double b, double macheps, double t, int
c=a;
fc=fa;
iter=1;
if (fabs(fc)<fabs(fb)){
if (std::abs(fc)<std::abs(fb)){
// Goto ext: from Brent ALGOL code
a=b;
b=c;
@@ -268,10 +268,10 @@ double Brent(FuncWrapper1D* f, double a, double b, double macheps, double t, int
d=b-a;
e=b-a;
m=0.5*(c-b);
tol=2*macheps*fabs(b)+t;
while (fabs(m)>tol && fb!=0){
tol=2*macheps*std::abs(b)+t;
while (std::abs(m)>tol && fb!=0){
// See if a bisection is forced
if (fabs(e)<tol || fabs(fa) <= fabs(fb)){
if (std::abs(e)<tol || std::abs(fa) <= std::abs(fb)){
m=0.5*(c-b);
d=e=m;
}
@@ -299,7 +299,7 @@ double Brent(FuncWrapper1D* f, double a, double b, double macheps, double t, int
s=e;
e=d;
m=0.5*(c-b);
if (2*p<3*m*q-fabs(tol*q) || p<fabs(0.5*s*q)){
if (2*p<3*m*q-std::abs(tol*q) || p<std::abs(0.5*s*q)){
d=p/q;
}
else{
@@ -309,7 +309,7 @@ double Brent(FuncWrapper1D* f, double a, double b, double macheps, double t, int
}
a=b;
fa=fb;
if (fabs(d)>tol){
if (std::abs(d)>tol){
b+=d;
}
else if (m>0){
@@ -328,7 +328,7 @@ double Brent(FuncWrapper1D* f, double a, double b, double macheps, double t, int
fc=fa;
d=e=b-a;
}
if (fabs(fc)<fabs(fb)){
if (std::abs(fc)<std::abs(fb)){
// Goto ext: from Brent ALGOL code
a=b;
b=c;
@@ -338,7 +338,7 @@ double Brent(FuncWrapper1D* f, double a, double b, double macheps, double t, int
fc=fa;
}
m=0.5*(c-b);
tol=2*macheps*fabs(b)+t;
tol=2*macheps*std::abs(b)+t;
iter+=1;
if (!ValidNumber(a)){
throw ValueError(format("Brent's method a is NAN").c_str());}
@@ -348,7 +348,7 @@ double Brent(FuncWrapper1D* f, double a, double b, double macheps, double t, int
throw ValueError(format("Brent's method c is NAN").c_str());}
if (iter>maxiter){
throw SolutionError(std::string("Brent's method reached maximum number of steps of %d ", maxiter));}
if (std::abs(fb)< 2*macheps*fabs(b)){
if (std::abs(fb)< 2*macheps*std::abs(b)){
return b;
}
}

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

@@ -284,7 +284,7 @@ TEST_CASE_METHOD(TransportValidationFixture, "Compare viscosities against publis
CHECK_NOTHROW(get_value(CoolProp::iviscosity));
CAPTURE(el.expected);
CAPTURE(actual);
CHECK(fabs(actual/el.expected-1) < el.tol);
CHECK(std::abs(actual/el.expected-1) < el.tol);
}
}
@@ -503,7 +503,7 @@ TEST_CASE_METHOD(TransportValidationFixture, "Compare thermal conductivities aga
get_value(CoolProp::iconductivity);
CAPTURE(el.expected);
CAPTURE(actual);
CHECK(fabs(actual/el.expected-1) < el.tol);
CHECK(std::abs(actual/el.expected-1) < el.tol);
}
}
@@ -606,8 +606,8 @@ public:
State.update(pair, x1, x2);
// Make sure we end up back at the same temperature and pressure we started out with
if(fabs(T-State.T()) > 1e-2) throw CoolProp::ValueError(format("Error on T [%g K] is greater than 1e-2",fabs(State.T()-T)));
if(fabs(p-State.p())/p*100 > 1e-2) throw CoolProp::ValueError(format("Error on p [%g %%] is greater than 1e-2 %%",fabs(p-State.p())/p ));
if(std::abs(T-State.T()) > 1e-2) throw CoolProp::ValueError(format("Error on T [%g K] is greater than 1e-2",std::abs(State.T()-T)));
if(std::abs(p-State.p())/p*100 > 1e-2) throw CoolProp::ValueError(format("Error on p [%g %%] is greater than 1e-2 %%",std::abs(p-State.p())/p ));
}
};
@@ -980,7 +980,7 @@ TEST_CASE("Triple point checks", "[triple_point]")
for (std::size_t i = 0; i < fluids.size(); ++i){
std::ostringstream ss1;
ss1 << "Triple point pressures matches for fluid " << fluids[i];
ss1 << "Triple point pressures matches for pure " << fluids[i];
SECTION(ss1.str(), "")
{
std::vector<std::string> names(1,fluids[i]);