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

Successive substitution is working again

Browse Source 
Signed-off-by: Ian Bell <ian.h.bell@gmail.com>
This commit is contained in:
Ian Bell
2014-08-25 15:37:30 +02:00
parent 9435b2d136
commit a6850b8f19
2 changed files with 27 additions and 8 deletions
Download Patch File Download Diff File Expand all files Collapse all files

31
src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp
View File

@@ -1590,16 +1590,15 @@ long double HelmholtzEOSMixtureBackend::solver_rho_Tp(long double T, long double
};
double call(double rhomolar){
this->rhomolar = rhomolar;
delta = rhomolar/rhor;
dalphar_dDelta = HEOS->calc_alphar_deriv_nocache(0, 1, HEOS->get_mole_fractions(), tau, delta);
peos = rhomolar*R_u*T*(1+delta*dalphar_dDelta);
delta = rhomolar/rhor; // needed for derivative
HEOS->update(DmolarT_INPUTS, rhomolar, T);
peos = HEOS->p();
r = (peos-p)/p;
return r;
};
double deriv(double rhomolar){
long double d2alphar_dDelta2 = HEOS->calc_alphar_deriv_nocache(0, 2, HEOS->get_mole_fractions(), tau, delta);
// dp/ddelta|tau / pspecified
return R_u*T*(1+2*delta*dalphar_dDelta+pow(delta, 2)*d2alphar_dDelta2)/p;
// dp/drho|T / pspecified
return R_u*T*(1+2*delta*HEOS->dalphar_dDelta()+pow(delta, 2)*HEOS->d2alphar_dDelta2())/p;
};
};
solver_TP_resid resid(this,T,p);
@@ -2032,7 +2031,25 @@ void HelmholtzEOSMixtureBackend::calc_all_alphar_deriv_cache(const std::vector<l
_d3alphar_dTau3 = derivs.d3alphar_dtau3;
}
else{
throw NotImplementedError("calc_all_alphar_deriv_nocache not yet implemented for mixtures");
std::size_t N = mole_fractions.size();
long double summer_base = 0, summer_dTau = 0, summer_dDelta = 0,
summer_dTau2 = 0, summer_dDelta2 = 0, summer_dDelta_dTau = 0;
for (std::size_t i = 0; i < N; ++i){
HelmholtzDerivatives derivs = components[i]->pEOS->alphar.all(tau, delta);
long double xi = mole_fractions[i];
summer_base += xi*derivs.alphar;
summer_dDelta += xi*derivs.dalphar_ddelta;
summer_dTau += xi*derivs.dalphar_dtau;
summer_dDelta2 += xi*derivs.d2alphar_ddelta2;
summer_dDelta_dTau += xi*derivs.d2alphar_ddelta_dtau;
summer_dTau2 += xi*derivs.d2alphar_dtau2;
}
_alphar = summer_base + Excess.alphar(tau, delta, mole_fractions);
_dalphar_dDelta = summer_dDelta + Excess.dalphar_dDelta(tau, delta, mole_fractions);
_dalphar_dTau = summer_dTau + Excess.dalphar_dTau(tau, delta, mole_fractions);
_d2alphar_dDelta2 = summer_dDelta2 + Excess.d2alphar_dDelta2(tau, delta, mole_fractions);
_d2alphar_dDelta_dTau = summer_dDelta_dTau + Excess.d2alphar_dDelta_dTau(tau, delta, mole_fractions);
_d2alphar_dTau2 = summer_dTau2 + Excess.d2alphar_dTau2(tau, delta, mole_fractions);
}
}

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

@@ -770,6 +770,8 @@ void SaturationSolvers::successive_substitution(HelmholtzEOSMixtureBackend &HEOS
}
x_and_y_from_K(beta, K, z, x, y);
normalize_vector(x);
normalize_vector(y);
HEOS.SatL->set_mole_fractions(x);
HEOS.SatV->set_mole_fractions(y);
@@ -779,7 +781,7 @@ void SaturationSolvers::successive_substitution(HelmholtzEOSMixtureBackend &HEOS
throw ValueError(format("saturation_p was unable to reach a solution within 50 iterations"));
}
}
while(std::abs(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);