Fixed(?) derivatives of spline function in two-phase derivs

This commit is contained in:
Ian Bell
2015-11-15 09:50:13 -07:00
parent c3e0e12c37
commit 9053983053
2 changed files with 24 additions and 6 deletions

View File

@@ -2966,16 +2966,17 @@ CoolPropDbl HelmholtzEOSMixtureBackend::calc_first_two_phase_deriv_splined(param
CoolPropDbl dhL_dp_sat = SatL->calc_first_saturation_deriv(h_key, p_key, *SatL, *SatV);
CoolPropDbl dhV_dp_sat = SatV->calc_first_saturation_deriv(h_key, p_key, *SatL, *SatV);
CoolPropDbl drhoL_dp_sat = SatL->calc_first_saturation_deriv(rho_key, p_key, *SatL, *SatV);
//CoolPropDbl drhoV_dp_sat = SatV->calc_first_saturation_deriv(rho_key, p_key, *SatL, *SatV);
//CoolPropDbl drho_dp_end = POW2(End->keyed_output(rho_key))*(x_end/POW2(rhoV)*drhoV_dp_sat + (1-x_end)/POW2(rhoL)*drhoL_dp_sat);
CoolPropDbl drhoV_dp_sat = SatV->calc_first_saturation_deriv(rho_key, p_key, *SatL, *SatV);
CoolPropDbl rhoV = SatV->keyed_output(rho_key);
CoolPropDbl rhoL = SatL->keyed_output(rho_key);
CoolPropDbl drho_dp_end = POW2(End->keyed_output(rho_key))*(x_end/POW2(rhoV)*drhoV_dp_sat + (1-x_end)/POW2(rhoL)*drhoL_dp_sat);
// Faking single-phase
//CoolPropDbl drho_dp__consth_liq = Liq->first_partial_deriv(rho_key, p_key, h_key);
CoolPropDbl d2rhodhdp_liq = Liq->second_partial_deriv(rho_key, h_key, p_key, p_key, h_key); // ?
// Derivatives at the end point
CoolPropDbl drho_dp__consth_end = End->calc_first_two_phase_deriv(rho_key, p_key, h_key);
// CoolPropDbl drho_dp__consth_end = End->calc_first_two_phase_deriv(rho_key, p_key, h_key);
CoolPropDbl d2rhodhdp_end = End->calc_second_two_phase_deriv(rho_key, h_key, p_key, p_key, h_key);
// Reminder:
@@ -2986,10 +2987,10 @@ CoolPropDbl HelmholtzEOSMixtureBackend::calc_first_two_phase_deriv_splined(param
// First pressure derivative at constant h of the coefficients a,b,c,d
// CoolPropDbl Abracket = (2*rho_liq - 2*rho_end + Delta_end * (drho_dh_liq__constp + drho_dh_end));
CoolPropDbl d_Abracket_dp_consth = (2*drhoL_dp_sat - 2*drho_dp__consth_end + Delta_end*(d2rhodhdp_liq + d2rhodhdp_end) + d_Delta_dp__consth*(drho_dh_liq__constp + drho_dh_end));
CoolPropDbl d_Abracket_dp_consth = (2*drhoL_dp_sat - 2*drho_dp_end + Delta_end*(d2rhodhdp_liq + d2rhodhdp_end) + d_Delta_end_dp__consth*(drho_dh_liq__constp + drho_dh_end));
CoolPropDbl da_dp = 1/POW3(Delta_end)*d_Abracket_dp_consth + Abracket*(-3/POW4(Delta_end)*d_Delta_end_dp__consth);
CoolPropDbl db_dp = - 6/POW3(Delta_end)*d_Delta_end_dp__consth*(rho_end - rho_liq)
+ (3/POW2(Delta_end))*(drho_dp__consth_end - drhoL_dp_sat)
+ 3/POW2(Delta_end)*(drho_dp_end - drhoL_dp_sat)
+ (1/POW2(Delta_end)*d_Delta_end_dp__consth) * (drho_dh_end + 2*drho_dh_liq__constp)
- (1/Delta_end) * (d2rhodhdp_end + 2*d2rhodhdp_liq);
CoolPropDbl dc_dp = d2rhodhdp_liq;

View File

@@ -1179,6 +1179,23 @@ TEST_CASE("Test second partial derivatives", "[derivatives]")
CHECK(std::abs((dudrho_T_num-dudrho_T_ana)/dudrho_T_ana) < tol);
CHECK(std::abs((d2udrho2_T_num-d2udrho2_T_ana)/d2udrho2_T_ana) < tol);
}
SECTION("Check second mixed partial(h,p) with respect to rho","")
{
shared_ptr<CoolProp::AbstractState> AS(CoolProp::AbstractState::factory("HEOS", "Propane"));
double dhmass = 1.0, T = 300;
AS->update(CoolProp::QT_INPUTS, 0.0, T);
double deriv1 = AS->first_partial_deriv(iDmass, iP, iHmass);
double deriv_analyt = AS->second_partial_deriv(iDmass, iP, iHmass, iHmass, iP);
double deriv_analyt2 = AS->second_partial_deriv(iDmass, iHmass, iP, iP, iHmass);
AS->update(CoolProp::HmassP_INPUTS, AS->hmass()-1, AS->p());
double deriv2 = AS->first_partial_deriv(iDmass, iP, iHmass);
double deriv_num = (deriv1-deriv2)/dhmass;
CAPTURE(deriv_num);
CAPTURE(deriv_analyt);
double tol = 1e-4;
CHECK(std::abs((deriv_num-deriv_analyt)/deriv_analyt) < tol);
}
}
TEST_CASE("REFPROP names for coolprop fluids", "[REFPROPName]")