diff --git a/src/Backends/Helmholtz/ExcessHEFunction.h b/src/Backends/Helmholtz/ExcessHEFunction.h index e45c16f0..2b7204f3 100644 --- a/src/Backends/Helmholtz/ExcessHEFunction.h +++ b/src/Backends/Helmholtz/ExcessHEFunction.h @@ -20,19 +20,26 @@ class DepartureFunction public: DepartureFunction(){}; virtual ~DepartureFunction(){}; + ResidualHelmholtzGeneralizedExponential phi; + HelmholtzDerivatives derivs; - /// The excess Helmholtz energy of the binary pair - /// Pure-virtual function (must be implemented in derived class - virtual double alphar(double tau, double delta) = 0; - virtual double dalphar_dDelta(double tau, double delta) = 0; - virtual double d2alphar_dDelta2(double tau, double delta) = 0; - virtual double d2alphar_dDelta_dTau(double tau, double delta) = 0; - virtual double dalphar_dTau(double tau, double delta) = 0; - virtual double d2alphar_dTau2(double tau, double delta) = 0; - virtual double d3alphar_dTau3(double tau, double delta) = 0; - virtual double d3alphar_dDelta_dTau2(double tau, double delta) = 0; - virtual double d3alphar_dDelta2_dTau(double tau, double delta) = 0; - virtual double d3alphar_dDelta3(double tau, double delta) = 0; + void update(double tau, double delta){ + derivs.reset(0.0); + phi.all(tau, delta, derivs); + }; + + double alphar(){ return derivs.alphar;}; + double dalphar_dDelta(){ return derivs.dalphar_ddelta;}; + double dalphar_dTau(){ return derivs.dalphar_dtau;}; + + double d2alphar_dDelta2(){return derivs.d2alphar_ddelta2;}; + double d2alphar_dDelta_dTau(){return derivs.d2alphar_ddelta_dtau;}; + double d2alphar_dTau2(){return derivs.d2alphar_dtau2;}; + + double d3alphar_dTau3(){ return derivs.d3alphar_dtau3; }; + double d3alphar_dDelta_dTau2(){ return derivs.d3alphar_ddelta_dtau2; }; + double d3alphar_dDelta2_dTau(){ return derivs.d3alphar_ddelta2_dtau; }; + double d3alphar_dDelta3(){ return derivs.d3alphar_ddelta3; }; }; /** \brief The departure function used by the GERG-2008 formulation @@ -44,10 +51,7 @@ public: * It is symmetric so \f$\alphar^r_{ij} = \alphar^r_{ji}\f$ */ class GERG2008DepartureFunction : public DepartureFunction -{ -protected: - bool using_gaussian; - ResidualHelmholtzGeneralizedExponential phi; +{ public: GERG2008DepartureFunction(){}; GERG2008DepartureFunction(const std::vector &n,const std::vector &d,const std::vector &t, @@ -64,11 +68,9 @@ public: } if (n.size() == Npower) { - using_gaussian = false; } else { - using_gaussian = true; std::vector _n(n.begin()+Npower, n.end()); std::vector _d(d.begin()+Npower, d.end()); std::vector _t(t.begin()+Npower, t.end()); @@ -80,17 +82,6 @@ public: } }; ~GERG2008DepartureFunction(){}; - - double alphar(double tau, double delta){return phi.base(tau, delta);}; - double dalphar_dDelta(double tau, double delta){return phi.dDelta(tau, delta);}; - double d2alphar_dDelta_dTau(double tau, double delta){return phi.dDelta_dTau(tau, delta);}; - double dalphar_dTau(double tau, double delta){return phi.dTau(tau, delta);}; - double d2alphar_dDelta2(double tau, double delta){return phi.dDelta2(tau, delta);}; - double d2alphar_dTau2(double tau, double delta){return phi.dTau2(tau, delta);}; - double d3alphar_dTau3(double tau, double delta){ return phi.dTau3(tau, delta); }; - double d3alphar_dDelta_dTau2(double tau, double delta){ return phi.dDelta_dTau2(tau, delta); }; - double d3alphar_dDelta2_dTau(double tau, double delta){ return phi.dDelta2_dTau(tau, delta); }; - double d3alphar_dDelta3(double tau, double delta){ return phi.dDelta3(tau, delta); }; }; /** \brief A polynomial/exponential departure function @@ -103,8 +94,6 @@ public: */ class ExponentialDepartureFunction : public DepartureFunction { -protected: - ResidualHelmholtzGeneralizedExponential phi; public: ExponentialDepartureFunction(){}; ExponentialDepartureFunction(const std::vector &n, const std::vector &d, @@ -117,19 +106,6 @@ public: phi.add_Power(_n, _d, _t, _l); }; ~ExponentialDepartureFunction(){}; - - double alphar(double tau, double delta){return phi.base(tau, delta);}; - double dalphar_dDelta(double tau, double delta){return phi.dDelta(tau, delta);}; - double d2alphar_dDelta_dTau(double tau, double delta){return phi.dDelta_dTau(tau, delta);}; - double dalphar_dTau(double tau, double delta){return phi.dTau(tau, delta);}; - double d2alphar_dDelta2(double tau, double delta){return phi.dDelta2(tau, delta);}; - double d2alphar_dTau2(double tau, double delta){return phi.dTau2(tau, delta);}; - - double d3alphar_dTau3(double tau, double delta){ return phi.dTau3(tau, delta); }; - double d3alphar_dDelta_dTau2(double tau, double delta){ return phi.dDelta_dTau2(tau, delta); }; - double d3alphar_dDelta2_dTau(double tau, double delta){ return phi.dDelta2_dTau(tau, delta); }; - double d3alphar_dDelta3(double tau, double delta){ return phi.dDelta3(tau, delta); }; - }; typedef shared_ptr DepartureFunctionPointer; @@ -152,233 +128,244 @@ public: DepartureFunctionMatrix[i].resize(N); } }; + /// Update the internal cached derivatives in each departure function + void update(double tau, double delta){ + for (std::size_t i = 0; i < N; i++){ + for (std::size_t j = i + 1; j < N; j++){ + DepartureFunctionMatrix[i][j]->update(tau, delta); + } + for (std::size_t j = 0; j < i; j++){ + DepartureFunctionMatrix[i][j]->update(tau, delta); + } + } + } - double alphar(double tau, double delta, const std::vector &x) + double alphar(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N-1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->alphar(tau,delta); + summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->alphar(); } } return summer; } - double dalphar_dDelta(double tau, double delta, const std::vector &x) + double dalphar_dDelta(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N-1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->dalphar_dDelta(tau,delta); + summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->dalphar_dDelta(); } } return summer; } - double d2alphar_dDelta2(double tau, double delta, const std::vector &x) + double d2alphar_dDelta2(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N-1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->d2alphar_dDelta2(tau,delta); + summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->d2alphar_dDelta2(); } } return summer; }; - double d2alphar_dDelta_dTau(double tau, double delta, const std::vector &x) + double d2alphar_dDelta_dTau(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N-1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->d2alphar_dDelta_dTau(tau,delta); + summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->d2alphar_dDelta_dTau(); } } return summer; } - double dalphar_dTau(double tau, double delta, const std::vector &x) + double dalphar_dTau(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N-1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->dalphar_dTau(tau,delta); + summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->dalphar_dTau(); } } return summer; }; - double d2alphar_dTau2(double tau, double delta, const std::vector &x) + double d2alphar_dTau2(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N-1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->d2alphar_dTau2(tau,delta); + summer += x[i]*x[j]*F[i][j]*DepartureFunctionMatrix[i][j]->d2alphar_dTau2(); } } return summer; }; - double d3alphar_dTau3(double tau, double delta, const std::vector &x) + double d3alphar_dTau3(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N - 1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i] * x[j] * F[i][j] * DepartureFunctionMatrix[i][j]->d3alphar_dTau3(tau, delta); + summer += x[i] * x[j] * F[i][j] * DepartureFunctionMatrix[i][j]->d3alphar_dTau3(); } } return summer; }; - double d3alphar_dDelta_dTau2(double tau, double delta, const std::vector &x) + double d3alphar_dDelta_dTau2(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N - 1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i] * x[j] * F[i][j] * DepartureFunctionMatrix[i][j]->d3alphar_dDelta_dTau2(tau, delta); + summer += x[i] * x[j] * F[i][j] * DepartureFunctionMatrix[i][j]->d3alphar_dDelta_dTau2(); } } return summer; }; - double d3alphar_dDelta2_dTau(double tau, double delta, const std::vector &x) + double d3alphar_dDelta2_dTau(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N - 1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i] * x[j] * F[i][j] * DepartureFunctionMatrix[i][j]->d3alphar_dDelta2_dTau(tau, delta); + summer += x[i] * x[j] * F[i][j] * DepartureFunctionMatrix[i][j]->d3alphar_dDelta2_dTau(); } } return summer; }; - double d3alphar_dDelta3(double tau, double delta, const std::vector &x) + double d3alphar_dDelta3(const std::vector &x) { double summer = 0; for (std::size_t i = 0; i < N - 1; i++) { for (std::size_t j = i + 1; j < N; j++) { - summer += x[i] * x[j] * F[i][j] * DepartureFunctionMatrix[i][j]->d3alphar_dDelta3(tau, delta); + summer += x[i] * x[j] * F[i][j] * DepartureFunctionMatrix[i][j]->d3alphar_dDelta3(); } } return summer; }; - double dalphar_dxi(double tau, double delta, const std::vector &x, std::size_t i) + double dalphar_dxi(const std::vector &x, std::size_t i) { double summer = 0; for (std::size_t k = 0; k < N; k++) { if (i != k) { - summer += x[k]*F[i][k]*DepartureFunctionMatrix[i][k]->alphar(tau,delta); + summer += x[k]*F[i][k]*DepartureFunctionMatrix[i][k]->alphar(); } } return summer; }; - double d2alphardxidxj(double tau, double delta, const std::vector &x, std::size_t i, std::size_t j) + double d2alphardxidxj(const std::vector &x, std::size_t i, std::size_t j) { if (i != j) { - return F[i][j]*DepartureFunctionMatrix[i][j]->alphar(tau,delta); + return F[i][j]*DepartureFunctionMatrix[i][j]->alphar(); } else { return 0; } }; - double d3alphar_dxi_dxj_dDelta(double tau, double delta, const std::vector &x, std::size_t i, std::size_t j) + double d3alphar_dxi_dxj_dDelta(const std::vector &x, std::size_t i, std::size_t j) { if (i != j) { - return F[i][j]*DepartureFunctionMatrix[i][j]->dalphar_dDelta(tau, delta); + return F[i][j]*DepartureFunctionMatrix[i][j]->dalphar_dDelta(); } else { return 0; } }; - double d3alphar_dxi_dxj_dTau(double tau, double delta, const std::vector &x, std::size_t i, std::size_t j) + double d3alphar_dxi_dxj_dTau(const std::vector &x, std::size_t i, std::size_t j) { if (i != j) { - return F[i][j]*DepartureFunctionMatrix[i][j]->dalphar_dTau(tau, delta); + return F[i][j]*DepartureFunctionMatrix[i][j]->dalphar_dTau(); } else { return 0; } }; - double d3alphardxidxjdxk(double tau, double delta, const std::vector &x, std::size_t i, std::size_t j, std::size_t k) + double d3alphardxidxjdxk(const std::vector &x, std::size_t i, std::size_t j, std::size_t k) { return 0; }; - double d2alphar_dxi_dTau(double tau, double delta, const std::vector &x, std::size_t i) + double d2alphar_dxi_dTau(const std::vector &x, std::size_t i) { double summer = 0; for (std::size_t k = 0; k < N; k++) { if (i != k) { - summer += x[k]*F[i][k]*DepartureFunctionMatrix[i][k]->dalphar_dTau(tau,delta); + summer += x[k]*F[i][k]*DepartureFunctionMatrix[i][k]->dalphar_dTau(); } } return summer; }; - double d2alphar_dxi_dDelta(double tau, double delta, const std::vector &x, std::size_t i) + double d2alphar_dxi_dDelta(const std::vector &x, std::size_t i) { double summer = 0; for (std::size_t k = 0; k < N; k++) { if (i != k) { - summer += x[k]*F[i][k]*DepartureFunctionMatrix[i][k]->dalphar_dDelta(tau,delta); + summer += x[k]*F[i][k]*DepartureFunctionMatrix[i][k]->dalphar_dDelta(); } } return summer; }; - double d3alphar_dxi_dDelta2(double tau, double delta, const std::vector &x, std::size_t i) + double d3alphar_dxi_dDelta2(const std::vector &x, std::size_t i) { double summer = 0; for (std::size_t k = 0; k < N; k++) { if (i != k) { - summer += x[k] * F[i][k] * DepartureFunctionMatrix[i][k]->d2alphar_dDelta2(tau, delta); + summer += x[k] * F[i][k] * DepartureFunctionMatrix[i][k]->d2alphar_dDelta2(); } } return summer; }; - double d3alphar_dxi_dTau2(double tau, double delta, const std::vector &x, std::size_t i) + double d3alphar_dxi_dTau2(const std::vector &x, std::size_t i) { double summer = 0; for (std::size_t k = 0; k < N; k++) { if (i != k) { - summer += x[k] * F[i][k] * DepartureFunctionMatrix[i][k]->d2alphar_dTau2(tau, delta); + summer += x[k] * F[i][k] * DepartureFunctionMatrix[i][k]->d2alphar_dTau2(); } } return summer; }; - double d3alphar_dxi_dDelta_dTau(double tau, double delta, const std::vector &x, std::size_t i) + double d3alphar_dxi_dDelta_dTau(const std::vector &x, std::size_t i) { double summer = 0; for (std::size_t k = 0; k < N; k++) { if (i != k) { - summer += x[k] * F[i][k] * DepartureFunctionMatrix[i][k]->d2alphar_dDelta_dTau(tau, delta); + summer += x[k] * F[i][k] * DepartureFunctionMatrix[i][k]->d2alphar_dDelta_dTau(); } } return summer; diff --git a/src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp b/src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp index b5168334..baabae08 100644 --- a/src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp +++ b/src/Backends/Helmholtz/HelmholtzEOSMixtureBackend.cpp @@ -943,7 +943,7 @@ void HelmholtzEOSMixtureBackend::pre_update(CoolProp::input_pairs &input_pair, C void HelmholtzEOSMixtureBackend::update(CoolProp::input_pairs input_pair, double value1, double value2 ) { if (get_debug_level() > 10){std::cout << format("%s (%d): update called with (%d: (%s), %g, %g)",__FILE__,__LINE__, input_pair, get_input_pair_short_desc(input_pair).c_str(), value1, value2) << std::endl;} - + CoolPropDbl ld_value1 = value1, ld_value2 = value2; pre_update(input_pair, ld_value1, ld_value2); value1 = ld_value1; value2 = ld_value2; @@ -2413,17 +2413,18 @@ void HelmholtzEOSMixtureBackend::calc_all_alphar_deriv_cache(const std::vector o(2); o[0] = L1; o[1] = M1; + double p = HEOS.p(); return o; }; std::vector > Jacobian(const std::vector &x) @@ -3016,13 +3018,13 @@ void HelmholtzEOSMixtureBackend::calc_critical_point(double rho0, double T0) J0(0,0) = (adjL*dLdTau).trace(); J0(0,1) = (adjL*dLdDelta).trace(); - //std::cout << J0 << std::endl; + std::cout << J0 << std::endl; std::vector r0 = call(x); // Build the Jacobian by column for (std::size_t i = 0; i < N; ++i) { xp = x; - epsilon = 0.0001; + epsilon = 0.00001; xp[i] += epsilon; r = call(xp); @@ -3036,6 +3038,7 @@ void HelmholtzEOSMixtureBackend::calc_critical_point(double rho0, double T0) }; }; Resid resid(*this); + CoolPropDbl Tr = T_reducing(); std::vector x, tau_delta(2); tau_delta[0] = T_reducing()/T0; tau_delta[1] = rho0/rhomolar_reducing(); std::vector td2 = tau_delta; std::string errstr; diff --git a/src/Backends/Helmholtz/MixtureDerivatives.cpp b/src/Backends/Helmholtz/MixtureDerivatives.cpp index c610c49a..6986d49c 100644 --- a/src/Backends/Helmholtz/MixtureDerivatives.cpp +++ b/src/Backends/Helmholtz/MixtureDerivatives.cpp @@ -5,19 +5,19 @@ namespace CoolProp{ CoolPropDbl MixtureDerivatives::dalphar_dxi(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return HEOS.components[i].EOS().baser(HEOS._tau, HEOS._delta) + HEOS.Excess.dalphar_dxi(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i); + return HEOS.components[i].EOS().baser(HEOS.tau(), HEOS.delta()) + HEOS.Excess.dalphar_dxi(HEOS.mole_fractions, i); } else if(xN_flag == XN_DEPENDENT){ std::vector &x = HEOS.mole_fractions; std::size_t N = x.size(); if (i == N-1) return 0; - double dar_dxi = HEOS.components[i].EOS().baser(HEOS._tau, HEOS._delta) - HEOS.components[N-1].EOS().baser(HEOS._tau, HEOS._delta); - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->alphar(HEOS._tau, HEOS._delta); + double dar_dxi = HEOS.components[i].EOS().baser(HEOS.tau(), HEOS.delta()) - HEOS.components[N-1].EOS().baser(HEOS.tau(), HEOS.delta()); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->alphar(); dar_dxi += (1-2*x[i])*FiNariN; for (std::size_t k = 0; k < N-1; ++k){ if (i == k) continue; - double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->alphar(HEOS._tau, HEOS._delta); - double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->alphar(HEOS._tau, HEOS._delta); + double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->alphar(); + double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->alphar(); dar_dxi += x[k]*(Fikarik - FiNariN - FkNarkN); } return dar_dxi; @@ -29,19 +29,19 @@ CoolPropDbl MixtureDerivatives::dalphar_dxi(HelmholtzEOSMixtureBackend &HEOS, st CoolPropDbl MixtureDerivatives::d2alphar_dxi_dTau(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return HEOS.components[i].EOS().dalphar_dTau(HEOS._tau, HEOS._delta) + HEOS.Excess.d2alphar_dxi_dTau(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i); + return HEOS.components[i].EOS().dalphar_dTau(HEOS._tau, HEOS._delta) + HEOS.Excess.d2alphar_dxi_dTau(HEOS.mole_fractions, i); } else if(xN_flag == XN_DEPENDENT){ std::vector &x = HEOS.mole_fractions; std::size_t N = x.size(); if (i==N-1) return 0; double d2ar_dxi_dTau = HEOS.components[i].EOS().dalphar_dTau(HEOS._tau, HEOS._delta) - HEOS.components[N-1].EOS().dalphar_dTau(HEOS._tau, HEOS._delta); - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->dalphar_dTau(HEOS._tau, HEOS._delta); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->dalphar_dTau(); d2ar_dxi_dTau += (1-2*x[i])*FiNariN; for (std::size_t k = 0; k < N-1; ++k){ if (i==k) continue; - double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->dalphar_dTau(HEOS._tau, HEOS._delta); - double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->dalphar_dTau(HEOS._tau, HEOS._delta); + double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->dalphar_dTau(); + double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->dalphar_dTau(); d2ar_dxi_dTau += x[k]*(Fikarik - FiNariN - FkNarkN); } return d2ar_dxi_dTau; @@ -54,19 +54,19 @@ CoolPropDbl MixtureDerivatives::d2alphar_dxi_dTau(HelmholtzEOSMixtureBackend &HE CoolPropDbl MixtureDerivatives::d2alphar_dxi_dDelta(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return HEOS.components[i].EOS().dalphar_dDelta(HEOS._tau, HEOS._delta) + HEOS.Excess.d2alphar_dxi_dDelta(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i); + return HEOS.components[i].EOS().dalphar_dDelta(HEOS.tau(), HEOS.delta()) + HEOS.Excess.d2alphar_dxi_dDelta(HEOS.mole_fractions, i); } else if(xN_flag == XN_DEPENDENT){ std::vector &x = HEOS.mole_fractions; std::size_t N = x.size(); if (i==N-1) return 0; - double d2ar_dxi_dDelta = HEOS.components[i].EOS().dalphar_dDelta(HEOS._tau, HEOS._delta) - HEOS.components[N-1].EOS().dalphar_dDelta(HEOS._tau, HEOS._delta); - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->dalphar_dDelta(HEOS._tau, HEOS._delta); + double d2ar_dxi_dDelta = HEOS.components[i].EOS().dalphar_dDelta(HEOS.tau(), HEOS.delta()) - HEOS.components[N-1].EOS().dalphar_dDelta(HEOS._tau, HEOS._delta); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->dalphar_dDelta(); d2ar_dxi_dDelta += (1-2*x[i])*FiNariN; for (std::size_t k = 0; k < N-1; ++k){ if (i==k) continue; - double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->dalphar_dDelta(HEOS._tau, HEOS._delta); - double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->dalphar_dDelta(HEOS._tau, HEOS._delta); + double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->dalphar_dDelta(); + double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->dalphar_dDelta(); d2ar_dxi_dDelta += x[k]*(Fikarik - FiNariN - FkNarkN); } return d2ar_dxi_dDelta; @@ -78,19 +78,19 @@ CoolPropDbl MixtureDerivatives::d2alphar_dxi_dDelta(HelmholtzEOSMixtureBackend & CoolPropDbl MixtureDerivatives::d3alphar_dxi_dDelta2(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return HEOS.components[i].EOS().d2alphar_dDelta2(HEOS._tau, HEOS._delta) + HEOS.Excess.d3alphar_dxi_dDelta2(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i); + return HEOS.components[i].EOS().d2alphar_dDelta2(HEOS.tau(), HEOS.delta()) + HEOS.Excess.d3alphar_dxi_dDelta2(HEOS.mole_fractions, i); } else if (xN_flag == XN_DEPENDENT){ std::vector &x = HEOS.mole_fractions; std::size_t N = x.size(); if (i==N-1) return 0; - double d3ar_dxi_dDelta2 = HEOS.components[i].EOS().d2alphar_dDelta2(HEOS._tau, HEOS._delta) - HEOS.components[N-1].EOS().d2alphar_dDelta2(HEOS._tau, HEOS._delta); - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->d2alphar_dDelta2(HEOS._tau, HEOS._delta); + double d3ar_dxi_dDelta2 = HEOS.components[i].EOS().d2alphar_dDelta2(HEOS.tau(), HEOS.delta()) - HEOS.components[N-1].EOS().d2alphar_dDelta2(HEOS.tau(), HEOS.delta()); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->d2alphar_dDelta2(); d3ar_dxi_dDelta2 += (1-2*x[i])*FiNariN; for (std::size_t k = 0; k < N-1; ++k){ if (i==k) continue; - double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->d2alphar_dDelta2(HEOS._tau, HEOS._delta); - double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->d2alphar_dDelta2(HEOS._tau, HEOS._delta); + double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->d2alphar_dDelta2(); + double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->d2alphar_dDelta2(); d3ar_dxi_dDelta2 += x[k]*(Fikarik - FiNariN - FkNarkN); } return d3ar_dxi_dDelta2; @@ -102,19 +102,19 @@ CoolPropDbl MixtureDerivatives::d3alphar_dxi_dDelta2(HelmholtzEOSMixtureBackend CoolPropDbl MixtureDerivatives::d3alphar_dxi_dTau2(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return HEOS.components[i].EOS().d2alphar_dTau2(HEOS._tau, HEOS._delta) + HEOS.Excess.d3alphar_dxi_dTau2(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i); + return HEOS.components[i].EOS().d2alphar_dTau2(HEOS.tau(), HEOS.delta()) + HEOS.Excess.d3alphar_dxi_dTau2(HEOS.mole_fractions, i); } else if (xN_flag == XN_DEPENDENT){ std::vector &x = HEOS.mole_fractions; std::size_t N = x.size(); if (i==N-1) return 0; - double d3ar_dxi_dTau2 = HEOS.components[i].EOS().d2alphar_dTau2(HEOS._tau, HEOS._delta) - HEOS.components[N-1].EOS().d2alphar_dTau2(HEOS._tau, HEOS._delta); - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->d2alphar_dTau2(HEOS._tau, HEOS._delta); + double d3ar_dxi_dTau2 = HEOS.components[i].EOS().d2alphar_dTau2(HEOS.tau(), HEOS.delta()) - HEOS.components[N-1].EOS().d2alphar_dTau2(HEOS.tau(), HEOS.delta()); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->d2alphar_dTau2(); d3ar_dxi_dTau2 += (1-2*x[i])*FiNariN; for (std::size_t k = 0; k < N-1; ++k){ if (i==k) continue; - double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->d2alphar_dTau2(HEOS._tau, HEOS._delta); - double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->d2alphar_dTau2(HEOS._tau, HEOS._delta); + double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->d2alphar_dTau2(); + double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->d2alphar_dTau2(); d3ar_dxi_dTau2 += x[k]*(Fikarik - FiNariN - FkNarkN); } return d3ar_dxi_dTau2; @@ -126,19 +126,19 @@ CoolPropDbl MixtureDerivatives::d3alphar_dxi_dTau2(HelmholtzEOSMixtureBackend &H CoolPropDbl MixtureDerivatives::d3alphar_dxi_dDelta_dTau(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return HEOS.components[i].EOS().d2alphar_dDelta_dTau(HEOS._tau, HEOS._delta) + HEOS.Excess.d3alphar_dxi_dDelta_dTau(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i); + return HEOS.components[i].EOS().d2alphar_dDelta_dTau(HEOS.tau(), HEOS.delta()) + HEOS.Excess.d3alphar_dxi_dDelta_dTau(HEOS.mole_fractions, i); } else if (xN_flag == XN_DEPENDENT){ std::vector &x = HEOS.mole_fractions; std::size_t N = x.size(); if (i==N-1) return 0; - double d3ar_dxi_dDelta_dTau = HEOS.components[i].EOS().d2alphar_dDelta_dTau(HEOS._tau, HEOS._delta) - HEOS.components[N-1].EOS().d2alphar_dDelta_dTau(HEOS._tau, HEOS._delta); - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->d2alphar_dDelta_dTau(HEOS._tau, HEOS._delta); + double d3ar_dxi_dDelta_dTau = HEOS.components[i].EOS().d2alphar_dDelta_dTau(HEOS.tau(), HEOS.delta()) - HEOS.components[N-1].EOS().d2alphar_dDelta_dTau(HEOS.tau(), HEOS.delta()); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->d2alphar_dDelta_dTau(); d3ar_dxi_dDelta_dTau += (1-2*x[i])*FiNariN; for (std::size_t k = 0; k < N-1; ++k){ if (i==k) continue; - double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->d2alphar_dDelta_dTau(HEOS._tau, HEOS._delta); - double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->d2alphar_dDelta_dTau(HEOS._tau, HEOS._delta); + double Fikarik = HEOS.Excess.F[i][k]*HEOS.Excess.DepartureFunctionMatrix[i][k]->d2alphar_dDelta_dTau(); + double FkNarkN = HEOS.Excess.F[k][N-1]*HEOS.Excess.DepartureFunctionMatrix[k][N-1]->d2alphar_dDelta_dTau(); d3ar_dxi_dDelta_dTau += x[k]*(Fikarik - FiNariN - FkNarkN); } return d3ar_dxi_dDelta_dTau; @@ -151,16 +151,16 @@ CoolPropDbl MixtureDerivatives::d3alphar_dxi_dDelta_dTau(HelmholtzEOSMixtureBack CoolPropDbl MixtureDerivatives::d2alphardxidxj(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return 0 + HEOS.Excess.d2alphardxidxj(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i, j); + return 0 + HEOS.Excess.d2alphardxidxj(HEOS.mole_fractions, i, j); } else if(xN_flag == XN_DEPENDENT){ std::size_t N = HEOS.mole_fractions.size(); if (i == N-1 || j == N-1){ return 0;} - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->alphar(HEOS._tau, HEOS._delta); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->alphar(); if (i == j) { return -2*FiNariN; } - double Fijarij = HEOS.Excess.F[i][j]*HEOS.Excess.DepartureFunctionMatrix[i][j]->alphar(HEOS._tau, HEOS._delta); - double FjNarjN = HEOS.Excess.F[j][N-1]*HEOS.Excess.DepartureFunctionMatrix[j][N-1]->alphar(HEOS._tau, HEOS._delta); + double Fijarij = HEOS.Excess.F[i][j]*HEOS.Excess.DepartureFunctionMatrix[i][j]->alphar(); + double FjNarjN = HEOS.Excess.F[j][N-1]*HEOS.Excess.DepartureFunctionMatrix[j][N-1]->alphar(); return Fijarij - FiNariN - FjNarjN; } else{ @@ -170,15 +170,15 @@ CoolPropDbl MixtureDerivatives::d2alphardxidxj(HelmholtzEOSMixtureBackend &HEOS, CoolPropDbl MixtureDerivatives::d3alphar_dxi_dxj_dDelta(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return 0 + HEOS.Excess.d3alphar_dxi_dxj_dDelta(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i, j); + return 0 + HEOS.Excess.d3alphar_dxi_dxj_dDelta(HEOS.mole_fractions, i, j); } else if (xN_flag == XN_DEPENDENT){ std::size_t N = HEOS.mole_fractions.size(); if (i == N-1 || j == N-1){ return 0; } - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->dalphar_dDelta(HEOS._tau, HEOS._delta); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->dalphar_dDelta(); if (i == j) { return -2*FiNariN; } - double Fijarij = HEOS.Excess.F[i][j]*HEOS.Excess.DepartureFunctionMatrix[i][j]->dalphar_dDelta(HEOS._tau, HEOS._delta); - double FjNarjN = HEOS.Excess.F[j][N-1]*HEOS.Excess.DepartureFunctionMatrix[j][N-1]->dalphar_dDelta(HEOS._tau, HEOS._delta); + double Fijarij = HEOS.Excess.F[i][j]*HEOS.Excess.DepartureFunctionMatrix[i][j]->dalphar_dDelta(); + double FjNarjN = HEOS.Excess.F[j][N-1]*HEOS.Excess.DepartureFunctionMatrix[j][N-1]->dalphar_dDelta(); return Fijarij - FiNariN - FjNarjN; } else{ @@ -188,15 +188,15 @@ CoolPropDbl MixtureDerivatives::d3alphar_dxi_dxj_dDelta(HelmholtzEOSMixtureBacke CoolPropDbl MixtureDerivatives::d3alphar_dxi_dxj_dTau(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return 0 + HEOS.Excess.d3alphar_dxi_dxj_dTau(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i, j); + return 0 + HEOS.Excess.d3alphar_dxi_dxj_dTau(HEOS.mole_fractions, i, j); } else if (xN_flag == XN_DEPENDENT){ std::size_t N = HEOS.mole_fractions.size(); if (i == N-1 || j == N-1){ return 0; } - double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->dalphar_dTau(HEOS._tau, HEOS._delta); + double FiNariN = HEOS.Excess.F[i][N-1]*HEOS.Excess.DepartureFunctionMatrix[i][N-1]->dalphar_dTau(); if (i == j) { return -2*FiNariN; } - double Fijarij = HEOS.Excess.F[i][j]*HEOS.Excess.DepartureFunctionMatrix[i][j]->dalphar_dTau(HEOS._tau, HEOS._delta); - double FjNarjN = HEOS.Excess.F[j][N-1]*HEOS.Excess.DepartureFunctionMatrix[j][N-1]->dalphar_dTau(HEOS._tau, HEOS._delta); + double Fijarij = HEOS.Excess.F[i][j]*HEOS.Excess.DepartureFunctionMatrix[i][j]->dalphar_dTau(); + double FjNarjN = HEOS.Excess.F[j][N-1]*HEOS.Excess.DepartureFunctionMatrix[j][N-1]->dalphar_dTau(); return Fijarij - FiNariN - FjNarjN; } else{ @@ -206,7 +206,7 @@ CoolPropDbl MixtureDerivatives::d3alphar_dxi_dxj_dTau(HelmholtzEOSMixtureBackend CoolPropDbl MixtureDerivatives::d3alphardxidxjdxk(HelmholtzEOSMixtureBackend &HEOS, std::size_t i, std::size_t j, std::size_t k, x_N_dependency_flag xN_flag) { if (xN_flag == XN_INDEPENDENT){ - return 0 + HEOS.Excess.d3alphardxidxjdxk(HEOS._tau, HEOS._delta, HEOS.mole_fractions, i, j, k); + return 0 + HEOS.Excess.d3alphardxidxjdxk(HEOS.mole_fractions, i, j, k); } else if (xN_flag == XN_DEPENDENT){ return 0;