diff --git a/include/CoolPropTools.h b/include/CoolPropTools.h index 055e2407..c9e21391 100644 --- a/include/CoolPropTools.h +++ b/include/CoolPropTools.h @@ -202,7 +202,7 @@ inline double min3(double x1, double x2, double x3){return std::min(std::min(x1, x2), x3);}; inline double max3(double x1, double x2, double x3){return std::max(std::max(x1, x2), x3);}; - inline bool double_equal(double a, double b){return fabs(a - b) <= 1 * DBL_EPSILON * std::max(fabs(a), fabs(b));}; + inline bool double_equal(double a, double b){return std::abs(a - b) <= 1 * DBL_EPSILON * std::max(std::abs(a), std::abs(b));}; template T max_abs_value(std::vector x) { @@ -210,7 +210,7 @@ std::size_t N = x.size(); for (std::size_t i = 0; i < N; ++i) { - T axi = fabs(x[i]); + T axi = std::abs(x[i]); if (axi > max){ max = axi; } } return max; @@ -263,11 +263,11 @@ /// Some functions related to testing and comparison of values bool inline check_abs(double A, double B, double D){ - double max = fabs(A); - double min = fabs(B); + double max = std::abs(A); + double min = std::abs(B); if (min>max) { max = min; - min = fabs(A); + min = std::abs(A); } if (max>DBL_EPSILON*1e3) return ( ( 1.0-min/max*1e0 ) < D ); else throw CoolProp::ValueError(format("Too small numbers: %f cannot be tested with an accepted error of %f for a machine precision of %f. ",max,D,DBL_EPSILON)); diff --git a/include/MatrixMath.h b/include/MatrixMath.h index cd37244a..382b8c0a 100644 --- a/include/MatrixMath.h +++ b/include/MatrixMath.h @@ -562,9 +562,9 @@ template size_t get_pivot_row(std::vector > *A, size_ for (size_t row = col; row < (*A).size(); row++) { val = (*A)[row][col]; - if (fabs(val) > max) + if (std::abs(val) > max) { - max = fabs(val); + max = std::abs(val); index = row; } } @@ -602,7 +602,7 @@ template std::vector > linsolve_Gauss_Jordan(std::vec // Find the pivot value pivot_row = get_pivot_row(&AB, col); - if (fabs(AB[pivot_row][col]) < 10*DBL_EPSILON){ throw ValueError(format("Zero occurred in row %d, the matrix is singular. ",pivot_row));} + if (std::abs(AB[pivot_row][col]) < 10*DBL_EPSILON){ throw ValueError(format("Zero occurred in row %d, the matrix is singular. ",pivot_row));} if (pivot_row>=col){ // Swap pivot row and current row @@ -671,7 +671,7 @@ template std::vector > linsolve_Gauss_Jordan(std::vec // pivot_row = 0; // pivot_element = 0.0; // for (size_t row = col; row < NrowA; row++){ -// tmp_element = fabs(AB[row][col]); +// tmp_element = std::abs(AB[row][col]); // if (tmp_element>pivot_element) { // pivot_element = tmp_element; // pivot_row = row; diff --git a/src/Backends/Helmholtz/FlashRoutines.cpp b/src/Backends/Helmholtz/FlashRoutines.cpp index f3758b48..7d7666e7 100644 --- a/src/Backends/Helmholtz/FlashRoutines.cpp +++ b/src/Backends/Helmholtz/FlashRoutines.cpp @@ -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)) { diff --git a/src/Backends/Helmholtz/Fluids/FluidLibrary.h b/src/Backends/Helmholtz/Fluids/FluidLibrary.h index adf457f0..d6343e1f 100644 --- a/src/Backends/Helmholtz/Fluids/FluidLibrary.h +++ b/src/Backends/Helmholtz/Fluids/FluidLibrary.h @@ -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 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 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); diff --git a/src/Backends/Helmholtz/TransportRoutines.cpp b/src/Backends/Helmholtz/TransportRoutines.cpp index 7aaa54d9..1ed7f941 100644 --- a/src/Backends/Helmholtz/TransportRoutines.cpp +++ b/src/Backends/Helmholtz/TransportRoutines.cpp @@ -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); diff --git a/src/Backends/Helmholtz/VLERoutines.cpp b/src/Backends/Helmholtz/VLERoutines.cpp index 67578a07..47dfb201 100644 --- a/src/Backends/Helmholtz/VLERoutines.cpp +++ b/src/Backends/Helmholtz/VLERoutines.cpp @@ -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");*/ diff --git a/src/Backends/REFPROP/REFPROPMixtureBackend.cpp b/src/Backends/REFPROP/REFPROPMixtureBackend.cpp index 1746f69c..94e3a45e 100644 --- a/src/Backends/REFPROP/REFPROPMixtureBackend.cpp +++ b/src/Backends/REFPROP/REFPROPMixtureBackend.cpp @@ -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); } } } diff --git a/src/CoolPropTools.cpp b/src/CoolPropTools.cpp index 44db51a6..020e7cc1 100644 --- a/src/CoolPropTools.cpp +++ b/src/CoolPropTools.cpp @@ -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 { diff --git a/src/Helmholtz.cpp b/src/Helmholtz.cpp index a5dc382d..396a7188 100644 --- a/src/Helmholtz.cpp +++ b/src/Helmholtz.cpp @@ -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); } } }; diff --git a/src/HumidAirProp.cpp b/src/HumidAirProp.cpp index fac7ceeb..134cdd7c 100644 --- a/src/HumidAirProp.cpp +++ b/src/HumidAirProp.cpp @@ -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); } } diff --git a/src/MatrixMath.cpp b/src/MatrixMath.cpp index e462093d..e2edc5e1 100644 --- a/src/MatrixMath.cpp +++ b/src/MatrixMath.cpp @@ -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 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 ); } } } diff --git a/src/PolyMath.cpp b/src/PolyMath.cpp index 3166ef08..3358b32c 100644 --- a/src/PolyMath.cpp +++ b/src/PolyMath.cpp @@ -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) NDNewtonRaphson_Jacobian(FuncWrapperND *f, std::vector f0,v,negative_f0; std::vector > 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)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){ + 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)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; } } diff --git a/src/Tests/CoolProp-Tests.cpp b/src/Tests/CoolProp-Tests.cpp index a0f17f16..c7537f3b 100644 --- a/src/Tests/CoolProp-Tests.cpp +++ b/src/Tests/CoolProp-Tests.cpp @@ -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 names(1,fluids[i]);