diff --git a/src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp b/src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp index d94c6200..7365cf33 100644 --- a/src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp +++ b/src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp @@ -1313,6 +1313,8 @@ void get_dT_drho(HelmholtzEOSMixtureBackend *HEOS, parameters index, long double dT = 1; drho = 0; break; case iDmolar: dT = 0; drho = 1; break; + case iDmass: + dT = 0; drho = HEOS->molar_mass(); break; case iP: { // dp/drho|T @@ -1321,30 +1323,51 @@ void get_dT_drho(HelmholtzEOSMixtureBackend *HEOS, parameters index, long double dT = rho*R*(1+delta*HEOS->dalphar_dDelta() - tau*delta*HEOS->d2alphar_dDelta_dTau()); break; } + case iHmass: case iHmolar: + { // dh/dT|rho dT = R*(-pow(tau,2)*(HEOS->d2alpha0_dTau2()+HEOS->d2alphar_dTau2()) + (1+delta*HEOS->dalphar_dDelta()-tau*delta*HEOS->d2alphar_dDelta_dTau())); - // dh/drho|T + // dh/drhomolar|T drho = T*R/rho*(tau*delta*HEOS->d2alphar_dDelta_dTau()+delta*HEOS->dalphar_dDelta()+pow(delta,2)*HEOS->d2alphar_dDelta2()); + if (index == iHmass){ + // dh/drho|T * drhomass/drhomolar where drhomass/drhomolar = mole mass + drho *= HEOS->molar_mass(); + } break; + } + case iSmass: case iSmolar: + { // ds/dT|rho dT = R/T*(-pow(tau,2)*(HEOS->d2alpha0_dTau2()+HEOS->d2alphar_dTau2())); // ds/drho|T drho = R/rho*(-(1+delta*HEOS->dalphar_dDelta()-tau*delta*HEOS->d2alphar_dDelta_dTau())); + if (index == iSmass){ + // ds/drho|T * drhomass/drhomolar where drhomass/drhomolar = mole mass + drho *= HEOS->molar_mass(); + } break; + } + case iUmass: case iUmolar: + { // du/dT|rho dT = R*(-pow(tau,2)*(HEOS->d2alpha0_dTau2()+HEOS->d2alphar_dTau2())); // du/drho|T drho = HEOS->T()*R/rho*(tau*delta*HEOS->d2alphar_dDelta_dTau()); + if (index == iUmass){ + // du/drho|T * drhomass/drhomolar where drhomass/drhomolar = mole mass + drho *= HEOS->molar_mass(); + } break; + } case iTau: dT = 1/dT_dtau; drho = 0; break; case iDelta: dT = 0; drho = 1/rhor; break; default: - throw ValueError(format("input to get_dtau_ddelta[%s] is invalid",get_parameter_information(index,"short").c_str())); + throw ValueError(format("input to get_dT_drho[%s] is invalid",get_parameter_information(index,"short").c_str())); } } void get_dT_drho_second_derivatives(HelmholtzEOSMixtureBackend *HEOS, int index, long double &dT2, long double &drho_dT, long double &drho2) @@ -1385,6 +1408,7 @@ void get_dT_drho_second_derivatives(HelmholtzEOSMixtureBackend *HEOS, int index, drho_dT = R*(1+2*delta*derivs.dalphar_ddelta +pow(delta,2)*derivs.d2alphar_ddelta2 - 2*delta*tau*derivs.d2alphar_ddelta_dtau - tau*pow(delta, 2)*derivs.d3alphar_ddelta2_dtau); break; } + case iHmass: case iHmolar: { // d2h/drho2|T @@ -1393,8 +1417,13 @@ void get_dT_drho_second_derivatives(HelmholtzEOSMixtureBackend *HEOS, int index, dT2 = R/T*pow(tau, 2)*(tau*(HEOS->d3alpha0_dTau3()+derivs.d3alphar_dtau3) + 2*(HEOS->d2alpha0_dTau2()+derivs.d2alphar_dtau2) + delta*derivs.d3alphar_ddelta_dtau2); // d2h/drho/dT drho_dT = R/rho*delta*(delta*derivs.d2alphar_ddelta2 - pow(tau,2)*derivs.d3alphar_ddelta_dtau2 + derivs.dalphar_ddelta - tau*delta*derivs.d3alphar_ddelta2_dtau - tau*derivs.d2alphar_ddelta_dtau); + if (index == iHmass){ + drho2 *= POW2(HEOS->molar_mass()); + drho_dT *= HEOS->molar_mass(); + } break; } + case iSmass: case iSmolar: { // d2s/rho2|T @@ -1403,8 +1432,13 @@ void get_dT_drho_second_derivatives(HelmholtzEOSMixtureBackend *HEOS, int index, dT2 = R*pow(tau/T, 2)*(tau*(HEOS->d3alpha0_dTau3()+derivs.d3alphar_dtau3)+3*(HEOS->d2alpha0_dTau2()+derivs.d2alphar_dtau2)); // d2s/drho/dT drho_dT = R/(T*rho)*(-pow(tau,2)*delta*derivs.d3alphar_ddelta_dtau2); + if (index == iSmass){ + drho2 *= POW2(HEOS->molar_mass()); + drho_dT *= HEOS->molar_mass(); + } break; } + case iUmass: case iUmolar: { // d2u/rho2|T @@ -1413,6 +1447,10 @@ void get_dT_drho_second_derivatives(HelmholtzEOSMixtureBackend *HEOS, int index, dT2 = R/T*pow(tau, 2)*(tau*(HEOS->d3alpha0_dTau3()+derivs.d3alphar_dtau3)+2*(HEOS->d2alpha0_dTau2()+derivs.d2alphar_dtau2)); // d2u/drho/dT drho_dT = R/rho*(-pow(tau,2)*delta*derivs.d3alphar_ddelta_dtau2); + if (index == iUmass){ + drho2 *= POW2(HEOS->molar_mass()); + drho_dT *= HEOS->molar_mass(); + } break; } default: diff --git a/src/Tests/CoolProp-Tests.cpp b/src/Tests/CoolProp-Tests.cpp index 6ae1a618..af6106d5 100644 --- a/src/Tests/CoolProp-Tests.cpp +++ b/src/Tests/CoolProp-Tests.cpp @@ -972,6 +972,24 @@ TEST_CASE("Test first partial derivatives using PropsSI", "[derivatives]") CHECK(rel_err_exact < 1e-7); CHECK(rel_err_approx < 1e-7); } + SECTION("Check dpdrho|T 3 ways for water using mass based","") + { + T = 80+273.15; + shared_ptr AS(CoolProp::AbstractState::factory("HEOS", "Water")); + AS->update(CoolProp::PT_INPUTS, 101325, T); + long double rhomass = AS->rhomass(); + double dpdrhomass__T_AbstractState = AS->first_partial_deriv(CoolProp::iP, CoolProp::iDmass, CoolProp::iT); + double dpdrhomass__T_PropsSI_num = (PropsSI("P","T",T,"Dmass",rhomass+1e-3,"Water") - PropsSI("P","T",T,"Dmass",rhomass-1e-3,"Water"))/(2*1e-3); + double dpdrhomass__T_PropsSI = PropsSI("d(P)/d(Dmass)|T","T",T,"P",101325,"Water"); + CAPTURE(rhomass); + CAPTURE(dpdrhomass__T_AbstractState); + CAPTURE(dpdrhomass__T_PropsSI_num); + CAPTURE(dpdrhomass__T_PropsSI); + double rel_err_exact = std::abs((dpdrhomass__T_AbstractState-dpdrhomass__T_PropsSI)/dpdrhomass__T_PropsSI); + double rel_err_approx = std::abs((dpdrhomass__T_PropsSI_num-dpdrhomass__T_PropsSI)/dpdrhomass__T_PropsSI); + CHECK(rel_err_exact < 1e-7); + CHECK(rel_err_approx < 1e-7); + } SECTION("Invalid first partial derivatives","") { CHECK(!ValidNumber(PropsSI("d()/d(P)|T","T",300,"P",101325,"n-Propane")));