CoolProp 7.0.0
An open-source fluid property and humid air property database
MixtureDerivatives.cpp
Go to the documentation of this file.
1#include "MixtureDerivatives.h"
2#include "Backends/Cubics/CubicBackend.h"
3
4namespace CoolProp {
5
6CoolPropDbl MixtureDerivatives::dln_fugacity_i_dT__constp_n(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
7 return dln_fugacity_coefficient_dT__constp_n(HEOS, i, xN_flag);
8}
9CoolPropDbl MixtureDerivatives::dln_fugacity_i_dp__constT_n(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
10 return dln_fugacity_coefficient_dp__constT_n(HEOS, i, xN_flag) + 1 / HEOS.p();
11}
12CoolPropDbl MixtureDerivatives::fugacity_i(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
13 return HEOS.mole_fractions[i] * HEOS.rhomolar() * HEOS.gas_constant() * HEOS.T() * exp(dnalphar_dni__constT_V_nj(HEOS, i, xN_flag));
14}
15CoolPropDbl MixtureDerivatives::ln_fugacity_coefficient(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
16 return HEOS.alphar() + ndalphar_dni__constT_V_nj(HEOS, i, xN_flag) - log(1 + HEOS._delta.pt() * HEOS.dalphar_dDelta());
17}
18CoolPropDbl MixtureDerivatives::dln_fugacity_i_dT__constrho_n(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
19 return 1 / HEOS.T() * (1 - HEOS.tau() * HEOS.dalphar_dTau() - HEOS.tau() * d_ndalphardni_dTau(HEOS, i, xN_flag));
20}
21CoolPropDbl MixtureDerivatives::dln_fugacity_i_drho__constT_n(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
22 return 1 / HEOS.rhomolar() * (1 + HEOS.delta() * HEOS.dalphar_dDelta() + HEOS.delta() * d_ndalphardni_dDelta(HEOS, i, xN_flag));
23}
24CoolPropDbl MixtureDerivatives::dnalphar_dni__constT_V_nj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
25 // GERG Equation 7.42
26 return HEOS.alphar() + ndalphar_dni__constT_V_nj(HEOS, i, xN_flag);
27}
28CoolPropDbl MixtureDerivatives::d2nalphar_dni_dT(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
29 return -HEOS._tau.pt() / HEOS._T * (HEOS.dalphar_dTau() + d_ndalphardni_dTau(HEOS, i, xN_flag));
30}
31CoolPropDbl MixtureDerivatives::dln_fugacity_coefficient_dT__constp_n(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
32 double T = HEOS._reducing.T / HEOS._tau.pt();
33 CoolPropDbl R_u = HEOS.gas_constant();
34 return d2nalphar_dni_dT(HEOS, i, xN_flag) + 1 / T - partial_molar_volume(HEOS, i, xN_flag) / (R_u * T) * dpdT__constV_n(HEOS);
35}
36CoolPropDbl MixtureDerivatives::partial_molar_volume(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
37 return -ndpdni__constT_V_nj(HEOS, i, xN_flag) / ndpdV__constT_n(HEOS);
38}
39
40CoolPropDbl MixtureDerivatives::dln_fugacity_coefficient_dp__constT_n(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
41 // GERG equation 7.30
42 CoolPropDbl R_u = HEOS.gas_constant();
43 double partial_molar_volumeval = partial_molar_volume(HEOS, i, xN_flag); // [m^3/mol]
44 double term1 = partial_molar_volumeval / (R_u * HEOS._T); // m^3/mol/(N*m)*mol = m^2/N = 1/Pa
45 double term2 = 1.0 / HEOS.p();
46 return term1 - term2;
47}
48CoolPropDbl MixtureDerivatives::dln_fugacity_i_dtau__constdelta_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
49 return -1 / HEOS.tau() + HEOS.dalphar_dTau() + d_ndalphardni_dTau(HEOS, i, xN_flag);
50}
51CoolPropDbl MixtureDerivatives::dln_fugacity_i_ddelta__consttau_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
52 return 1 + HEOS.delta() * HEOS.dalphar_dDelta() + HEOS.delta() * d_ndalphardni_dDelta(HEOS, i, xN_flag);
53}
54CoolPropDbl MixtureDerivatives::dln_fugacity_dxj__constT_p_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
55 x_N_dependency_flag xN_flag) {
56 // This is a term to which some more might be added depending on i and j
57 CoolPropDbl val = dln_fugacity_coefficient_dxj__constT_p_xi(HEOS, i, j, xN_flag);
58 const std::vector<CoolPropDbl>& x = HEOS.get_mole_fractions();
59 std::size_t N = x.size();
60 if (i == N - 1) {
61 val += -1 / x[N - 1];
62 } else if (i == j) {
63 val += 1 / x[j];
64 }
65 return val;
66}
67CoolPropDbl MixtureDerivatives::dln_fugacity_dxj__constT_rho_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
68 x_N_dependency_flag xN_flag) {
69 if (xN_flag == XN_INDEPENDENT) {
70 throw ValueError("dln_fugacity_dxj__constT_rho_xi only valid for xN_DEPENDENT for now");
71 }
72 CoolPropDbl rhor = HEOS.Reducing->rhormolar(HEOS.get_mole_fractions());
73 CoolPropDbl Tr = HEOS.Reducing->Tr(HEOS.get_mole_fractions());
74 CoolPropDbl dTrdxj = HEOS.Reducing->dTrdxi__constxj(HEOS.get_mole_fractions(), j, xN_flag);
75 CoolPropDbl drhordxj = HEOS.Reducing->drhormolardxi__constxj(HEOS.get_mole_fractions(), j, xN_flag);
76
77 // These lines are all the same
78 CoolPropDbl line1 = dln_fugacity_i_dtau__constdelta_x(HEOS, i, xN_flag) * 1 / HEOS.T() * dTrdxj;
79 CoolPropDbl line2 = -dln_fugacity_i_ddelta__consttau_x(HEOS, i, xN_flag) * 1 / rhor * drhordxj;
80 CoolPropDbl line4 = HEOS.residual_helmholtz->dalphar_dxi(HEOS, j, xN_flag) + d_ndalphardni_dxj__constdelta_tau_xi(HEOS, i, j, xN_flag);
81
82 const std::vector<CoolPropDbl>& x = HEOS.get_mole_fractions();
83 std::size_t N = x.size();
84
85 CoolPropDbl line3 = 1 / rhor * HEOS.Reducing->drhormolardxi__constxj(x, j, xN_flag) + 1 / Tr * HEOS.Reducing->dTrdxi__constxj(x, j, xN_flag);
86 ;
87
88 // This is a term to which some more might be added depending on i and j
89 if (i == N - 1) {
90 line3 += -1 / x[N - 1];
91 } else if (i == j) {
92 line3 += 1 / x[j];
93 } else {
94 }
95
96 return line1 + line2 + line3 + line4;
97}
98
99CoolPropDbl MixtureDerivatives::ndln_fugacity_i_dnj__constT_V_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
100 x_N_dependency_flag xN_flag) {
101 double s = (HEOS.mole_fractions[i] > DBL_EPSILON) ? Kronecker_delta(i, j) / HEOS.mole_fractions[i] : 0;
102 return s + nd2nalphardnidnj__constT_V(HEOS, i, j, xN_flag);
103}
104CoolPropDbl MixtureDerivatives::d_ndln_fugacity_i_dnj_dtau__constdelta_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
105 x_N_dependency_flag xN_flag) {
106 return d_ndalphardni_dTau(HEOS, j, xN_flag) + d_nd_ndalphardni_dnj_dTau__constdelta_x(HEOS, i, j, xN_flag);
107}
108CoolPropDbl MixtureDerivatives::d2_ndln_fugacity_i_dnj_dtau2__constdelta_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
109 x_N_dependency_flag xN_flag) {
110 return d2_ndalphardni_dTau2(HEOS, j, xN_flag) + d2_nd_ndalphardni_dnj_dTau2__constdelta_x(HEOS, i, j, xN_flag);
111}
112CoolPropDbl MixtureDerivatives::d2_ndln_fugacity_i_dnj_ddelta2__consttau_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
113 x_N_dependency_flag xN_flag) {
114 return d2_ndalphardni_dDelta2(HEOS, j, xN_flag) + d2_nd_ndalphardni_dnj_dDelta2__consttau_x(HEOS, i, j, xN_flag);
115}
116CoolPropDbl MixtureDerivatives::d2_ndln_fugacity_i_dnj_ddelta_dtau__constx(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
117 x_N_dependency_flag xN_flag) {
118 return d2_ndalphardni_dDelta_dTau(HEOS, j, xN_flag) + d2_nd_ndalphardni_dnj_dDelta_dTau__constx(HEOS, i, j, xN_flag);
119}
120
121CoolPropDbl MixtureDerivatives::d_ndln_fugacity_i_dnj_ddelta__consttau_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
122 x_N_dependency_flag xN_flag) {
123 return d_ndalphardni_dDelta(HEOS, j, xN_flag) + d_nd_ndalphardni_dnj_dDelta__consttau_x(HEOS, i, j, xN_flag);
124}
125CoolPropDbl MixtureDerivatives::d_ndln_fugacity_i_dnj_ddxk__consttau_delta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
126 std::size_t k, x_N_dependency_flag xN_flag) {
127 CoolPropDbl s = (HEOS.mole_fractions[i] > DBL_EPSILON) ? -Kronecker_delta(i, j) * Kronecker_delta(i, k) / pow(HEOS.mole_fractions[i], 2) : 0;
128 return s + d_ndalphardni_dxj__constdelta_tau_xi(HEOS, j, k, xN_flag) + d_nd_ndalphardni_dnj_dxk__consttau_delta(HEOS, i, j, k, xN_flag);
129}
130CoolPropDbl MixtureDerivatives::d2_ndln_fugacity_i_dnj_dxk_dTau__constdelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
131 std::size_t k, x_N_dependency_flag xN_flag) {
132 return d2_ndalphardni_dxj_dTau__constdelta_xi(HEOS, j, k, xN_flag) + d2_nd_ndalphardni_dnj_dxk_dTau__constdelta(HEOS, i, j, k, xN_flag);
133}
134CoolPropDbl MixtureDerivatives::d2_ndln_fugacity_i_dnj_dxk_dDelta__consttau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
135 std::size_t k, x_N_dependency_flag xN_flag) {
136 return d2_ndalphardni_dxj_dDelta__consttau_xi(HEOS, j, k, xN_flag) + d2_nd_ndalphardni_dnj_dxk_dDelta__consttau(HEOS, i, j, k, xN_flag);
137}
138CoolPropDbl MixtureDerivatives::nd_ndln_fugacity_i_dnj_dnk__constT_V_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, std::size_t k,
139 x_N_dependency_flag xN_flag) {
140 double sum = d_ndln_fugacity_i_dnj_dtau__constdelta_x(HEOS, i, j, xN_flag) * ndtaudni__constT_V_nj(HEOS, k, xN_flag)
141 + d_ndln_fugacity_i_dnj_ddelta__consttau_x(HEOS, i, j, xN_flag) * nddeltadni__constT_V_nj(HEOS, k, xN_flag)
142 + d_ndln_fugacity_i_dnj_ddxk__consttau_delta(HEOS, i, j, k, xN_flag);
143 std::size_t mmax = HEOS.mole_fractions.size();
144 if (xN_flag == XN_DEPENDENT) {
145 mmax--;
146 }
147 for (unsigned int m = 0; m < mmax; ++m) {
148 sum -= HEOS.mole_fractions[m] * d_ndln_fugacity_i_dnj_ddxk__consttau_delta(HEOS, i, j, m, xN_flag);
149 }
150 return sum;
151}
152
153CoolPropDbl MixtureDerivatives::dln_fugacity_coefficient_dxj__constT_p_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
154 x_N_dependency_flag xN_flag) {
155 // Gernert 3.115
156 CoolPropDbl R_u = HEOS.gas_constant();
157 // partial molar volume is -dpdn/dpdV, so need to flip the sign here
158 return d2nalphar_dxj_dni__constT_V(HEOS, j, i, xN_flag)
159 - partial_molar_volume(HEOS, i, xN_flag) / (R_u * HEOS._T) * dpdxj__constT_V_xi(HEOS, j, xN_flag);
160}
161CoolPropDbl MixtureDerivatives::dpdxj__constT_V_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t j, x_N_dependency_flag xN_flag) {
162 // Gernert 3.130
163 CoolPropDbl R_u = HEOS.gas_constant();
164 return HEOS._rhomolar * R_u * HEOS._T
165 * (ddelta_dxj__constT_V_xi(HEOS, j, xN_flag) * HEOS.dalphar_dDelta()
166 + HEOS._delta.pt() * d_dalpharddelta_dxj__constT_V_xi(HEOS, j, xN_flag));
167}
168
169CoolPropDbl MixtureDerivatives::d_dalpharddelta_dxj__constT_V_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t j, x_N_dependency_flag xN_flag) {
170 // Gernert Equation 3.134 (Catch test provided)
171 return HEOS.d2alphar_dDelta2() * ddelta_dxj__constT_V_xi(HEOS, j, xN_flag) + HEOS.d2alphar_dDelta_dTau() * dtau_dxj__constT_V_xi(HEOS, j, xN_flag)
172 + HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, j, xN_flag);
173}
174
175CoolPropDbl MixtureDerivatives::dalphar_dxj__constT_V_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t j, x_N_dependency_flag xN_flag) {
176 //Gernert 3.119 (Catch test provided)
177 return HEOS.dalphar_dDelta() * ddelta_dxj__constT_V_xi(HEOS, j, xN_flag) + HEOS.dalphar_dTau() * dtau_dxj__constT_V_xi(HEOS, j, xN_flag)
178 + HEOS.residual_helmholtz->dalphar_dxi(HEOS, j, xN_flag);
179}
180CoolPropDbl MixtureDerivatives::d_ndalphardni_dxj__constT_V_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
181 x_N_dependency_flag xN_flag) {
182 // Gernert 3.118
183 return d_ndalphardni_dxj__constdelta_tau_xi(HEOS, i, j, xN_flag)
184 + ddelta_dxj__constT_V_xi(HEOS, j, xN_flag) * d_ndalphardni_dDelta(HEOS, i, xN_flag)
185 + dtau_dxj__constT_V_xi(HEOS, j, xN_flag) * d_ndalphardni_dTau(HEOS, i, xN_flag);
186}
187CoolPropDbl MixtureDerivatives::ddelta_dxj__constT_V_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t j, x_N_dependency_flag xN_flag) {
188 // Gernert 3.121 (Catch test provided)
189 return -HEOS._delta.pt() / HEOS._reducing.rhomolar * HEOS.Reducing->drhormolardxi__constxj(HEOS.mole_fractions, j, xN_flag);
190}
191CoolPropDbl MixtureDerivatives::dtau_dxj__constT_V_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t j, x_N_dependency_flag xN_flag) {
192 // Gernert 3.122 (Catch test provided)
193 return 1 / HEOS._T * HEOS.Reducing->dTrdxi__constxj(HEOS.mole_fractions, j, xN_flag);
194}
195
196CoolPropDbl MixtureDerivatives::dpdT__constV_n(HelmholtzEOSMixtureBackend& HEOS) {
197 CoolPropDbl R_u = HEOS.gas_constant();
198 return HEOS._rhomolar * R_u * (1 + HEOS._delta.pt() * HEOS.dalphar_dDelta() - HEOS._delta.pt() * HEOS._tau.pt() * HEOS.d2alphar_dDelta_dTau());
199}
200CoolPropDbl MixtureDerivatives::dpdrho__constT_n(HelmholtzEOSMixtureBackend& HEOS) {
201 CoolPropDbl R_u = HEOS.gas_constant();
202 return R_u * HEOS._T * (1 + 2 * HEOS._delta.pt() * HEOS.dalphar_dDelta() + pow(HEOS._delta.pt(), 2) * HEOS.d2alphar_dDelta2());
203}
204CoolPropDbl MixtureDerivatives::ndpdV__constT_n(HelmholtzEOSMixtureBackend& HEOS) {
205 CoolPropDbl R_u = HEOS.gas_constant();
206 return -pow(HEOS._rhomolar, 2) * R_u * HEOS._T
207 * (1 + 2 * HEOS._delta.pt() * HEOS.dalphar_dDelta() + pow(HEOS._delta.pt(), 2) * HEOS.d2alphar_dDelta2());
208}
209CoolPropDbl MixtureDerivatives::ndpdni__constT_V_nj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
210 // Eqn 7.64 and 7.63
211 CoolPropDbl R_u = HEOS.gas_constant();
212 double ndrhorbar_dni__constnj = HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag);
213 double ndTr_dni__constnj = HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag);
214 double summer = 0;
215 std::size_t kmax = HEOS.mole_fractions.size();
216 if (xN_flag == XN_DEPENDENT) {
217 kmax--;
218 }
219 for (unsigned int k = 0; k < kmax; ++k) {
220 summer += HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, k, xN_flag);
221 }
222 double nd2alphar_dni_dDelta = HEOS._delta.pt() * HEOS.d2alphar_dDelta2() * (1 - 1 / HEOS._reducing.rhomolar * ndrhorbar_dni__constnj)
223 + HEOS._tau.pt() * HEOS.d2alphar_dDelta_dTau() / HEOS._reducing.T * ndTr_dni__constnj
224 + HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, i, xN_flag) - summer;
225 return HEOS._rhomolar * R_u * HEOS._T
226 * (1 + HEOS._delta.pt() * HEOS.dalphar_dDelta() * (2 - 1 / HEOS._reducing.rhomolar * ndrhorbar_dni__constnj)
227 + HEOS._delta.pt() * nd2alphar_dni_dDelta);
228}
229
230CoolPropDbl MixtureDerivatives::ndalphar_dni__constT_V_nj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
231 double term1 = HEOS._delta.pt() * HEOS.dalphar_dDelta()
232 * (1 - 1 / HEOS._reducing.rhomolar * HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag));
233 double term2 = HEOS._tau.pt() * HEOS.dalphar_dTau() * (1 / HEOS._reducing.T) * HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag);
234
235 double s = 0;
236 std::size_t kmax = HEOS.mole_fractions.size();
237 if (xN_flag == XN_DEPENDENT) {
238 kmax--;
239 }
240 for (unsigned int k = 0; k < kmax; k++) {
241 s += HEOS.mole_fractions[k] * HEOS.residual_helmholtz->dalphar_dxi(HEOS, k, xN_flag);
242 }
243 double term3 = HEOS.residual_helmholtz->dalphar_dxi(HEOS, i, xN_flag);
244 return term1 + term2 + term3 - s;
245}
246CoolPropDbl MixtureDerivatives::ndln_fugacity_coefficient_dnj__constT_p(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
247 x_N_dependency_flag xN_flag) {
248 CoolPropDbl R_u = HEOS.gas_constant();
249 return nd2nalphardnidnj__constT_V(HEOS, j, i, xN_flag) + 1
250 - partial_molar_volume(HEOS, j, xN_flag) / (R_u * HEOS._T) * ndpdni__constT_V_nj(HEOS, i, xN_flag);
251}
252CoolPropDbl MixtureDerivatives::nddeltadni__constT_V_nj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
253 return HEOS._delta.pt() - HEOS._delta.pt() / HEOS._reducing.rhomolar * HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag);
254}
255CoolPropDbl MixtureDerivatives::d_nddeltadni_dDelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
256 return 1 - 1 / HEOS._reducing.rhomolar * HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag);
257}
258CoolPropDbl MixtureDerivatives::d_nddeltadni_dxj__constdelta_tau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
259 x_N_dependency_flag xN_flag) {
260 double rhor = HEOS._reducing.rhomolar;
261 return -HEOS.delta() / rhor
262 * (HEOS.Reducing->d_ndrhorbardni_dxj__constxi(HEOS.mole_fractions, i, j, xN_flag)
263 - 1 / rhor * HEOS.Reducing->drhormolardxi__constxj(HEOS.mole_fractions, j, xN_flag)
264 * HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag));
265}
266CoolPropDbl MixtureDerivatives::d2_nddeltadni_dxj_dDelta__consttau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
267 x_N_dependency_flag xN_flag) {
268 return d_nddeltadni_dxj__constdelta_tau(HEOS, i, j, xN_flag) / HEOS.delta();
269}
270
271CoolPropDbl MixtureDerivatives::ndtaudni__constT_V_nj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
272 return HEOS._tau.pt() / HEOS._reducing.T * HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag);
273}
274
275CoolPropDbl MixtureDerivatives::d_ndtaudni_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
276 return 1 / HEOS._reducing.T * HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag);
277}
278CoolPropDbl MixtureDerivatives::d_ndtaudni_dxj__constdelta_tau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
279 x_N_dependency_flag xN_flag) {
280 double Tr = HEOS._reducing.T;
281 return HEOS.tau() / Tr
282 * (HEOS.Reducing->d_ndTrdni_dxj__constxi(HEOS.mole_fractions, i, j, xN_flag)
283 - 1 / Tr * HEOS.Reducing->dTrdxi__constxj(HEOS.mole_fractions, j, xN_flag)
284 * HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag));
285}
286CoolPropDbl MixtureDerivatives::d2_ndtaudni_dxj_dTau__constdelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
287 x_N_dependency_flag xN_flag) {
288 return d_ndtaudni_dxj__constdelta_tau(HEOS, i, j, xN_flag) / HEOS.tau();
289}
290CoolPropDbl MixtureDerivatives::d_ndalphardni_dxj__constdelta_tau_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
291 x_N_dependency_flag xN_flag) {
292 double line1 = HEOS._delta.pt() * HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, j, xN_flag)
293 * (1 - 1 / HEOS._reducing.rhomolar * HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag));
294 double line3 = HEOS._tau.pt() * HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, j, xN_flag) * (1 / HEOS._reducing.T)
295 * HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag);
296 double line2 = -HEOS._delta.pt() * HEOS.dalphar_dDelta() * (1 / HEOS._reducing.rhomolar)
297 * (HEOS.Reducing->d_ndrhorbardni_dxj__constxi(HEOS.mole_fractions, i, j, xN_flag)
298 - 1 / HEOS._reducing.rhomolar * HEOS.Reducing->drhormolardxi__constxj(HEOS.mole_fractions, j, xN_flag)
299 * HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag));
300 double line4 = HEOS._tau.pt() * HEOS.dalphar_dTau() * (1 / HEOS._reducing.T)
301 * (HEOS.Reducing->d_ndTrdni_dxj__constxi(HEOS.mole_fractions, i, j, xN_flag)
302 - 1 / HEOS._reducing.T * HEOS.Reducing->dTrdxi__constxj(HEOS.mole_fractions, j, xN_flag)
303 * HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag));
304
305 double s = 0;
306 std::size_t kmax = HEOS.mole_fractions.size();
307 if (xN_flag == XN_DEPENDENT) {
308 kmax--;
309 }
310 for (unsigned int k = 0; k < kmax; k++) {
311 s += HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d2alphardxidxj(HEOS, j, k, xN_flag);
312 }
313 double line5 = HEOS.residual_helmholtz->d2alphardxidxj(HEOS, i, j, xN_flag) - HEOS.residual_helmholtz->dalphar_dxi(HEOS, j, xN_flag) - s;
314 return line1 + line2 + line3 + line4 + line5;
315}
316
317CoolPropDbl MixtureDerivatives::nd2nalphardnidnj__constT_V(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
318 x_N_dependency_flag xN_flag) {
319 return ndalphar_dni__constT_V_nj(HEOS, j, xN_flag) // First term from 7.46
320 + nd_ndalphardni_dnj__constT_V(HEOS, i, j, xN_flag); // 7.47
321}
322
323CoolPropDbl MixtureDerivatives::nd_ndalphardni_dnj__constT_V(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
324 x_N_dependency_flag xN_flag) {
325 double line1 = d_ndalphardni_dDelta(HEOS, i, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag);
326 double line2 = d_ndalphardni_dTau(HEOS, i, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag);
327 double line3 = d_ndalphardni_dxj__constdelta_tau_xi(HEOS, i, j, xN_flag);
328 std::size_t kmax = HEOS.mole_fractions.size();
329 if (xN_flag == XN_DEPENDENT) {
330 kmax--;
331 }
332 for (unsigned int k = 0; k < kmax; k++) {
333 line3 -= HEOS.mole_fractions[k] * d_ndalphardni_dxj__constdelta_tau_xi(HEOS, i, k, xN_flag);
334 }
335 return line1 + line2 + line3;
336}
337CoolPropDbl MixtureDerivatives::d_nd_ndalphardni_dnj_dTau__constdelta_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
338 x_N_dependency_flag xN_flag) {
339 double line1 = d2_ndalphardni_dDelta_dTau(HEOS, i, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag);
340 double line2 = d2_ndalphardni_dTau2(HEOS, i, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag)
341 + d_ndalphardni_dTau(HEOS, i, xN_flag) * d_ndtaudni_dTau(HEOS, j, xN_flag);
342 double summer = 0;
343 std::size_t kmax = HEOS.mole_fractions.size();
344 if (xN_flag == XN_DEPENDENT) {
345 kmax--;
346 }
347 for (unsigned int k = 0; k < kmax; k++) {
348 summer += HEOS.mole_fractions[k] * d2_ndalphardni_dxj_dTau__constdelta_xi(HEOS, i, k, xN_flag);
349 }
350 double line3 = d2_ndalphardni_dxj_dTau__constdelta_xi(HEOS, i, j, xN_flag) - summer;
351 return line1 + line2 + line3;
352}
353CoolPropDbl MixtureDerivatives::d2_nd_ndalphardni_dnj_dTau2__constdelta_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
354 x_N_dependency_flag xN_flag) {
355 double line1 = d3_ndalphardni_dDelta_dTau2(HEOS, i, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag);
356 double line2 = 2 * d2_ndalphardni_dTau2(HEOS, i, xN_flag) * d_ndtaudni_dTau(HEOS, j, xN_flag);
357 double line3 = d3_ndalphardni_dTau3(HEOS, i, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag);
358 double summer = 0;
359 std::size_t kmax = HEOS.mole_fractions.size();
360 if (xN_flag == XN_DEPENDENT) {
361 kmax--;
362 }
363 for (unsigned int k = 0; k < kmax; k++) {
364 summer += HEOS.mole_fractions[k] * d3_ndalphardni_dxj_dTau2__constdelta_xi(HEOS, i, k, xN_flag);
365 }
366 double line4 = d3_ndalphardni_dxj_dTau2__constdelta_xi(HEOS, i, j, xN_flag) - summer;
367 return line1 + line2 + line3 + line4;
368}
369CoolPropDbl MixtureDerivatives::d2_nd_ndalphardni_dnj_dDelta2__consttau_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
370 x_N_dependency_flag xN_flag) {
371 double line1 = d3_ndalphardni_dDelta3(HEOS, i, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag);
372 double line2 = 2 * d2_ndalphardni_dDelta2(HEOS, i, xN_flag) * d_nddeltadni_dDelta(HEOS, j, xN_flag);
373 double line3 = d3_ndalphardni_dDelta2_dTau(HEOS, i, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag);
374 double summer = 0;
375 std::size_t kmax = HEOS.mole_fractions.size();
376 if (xN_flag == XN_DEPENDENT) {
377 kmax--;
378 }
379 for (unsigned int k = 0; k < kmax; k++) {
380 summer += HEOS.mole_fractions[k] * d3_ndalphardni_dxj_dDelta2__consttau_xi(HEOS, i, k, xN_flag);
381 }
382 double line4 = d3_ndalphardni_dxj_dDelta2__consttau_xi(HEOS, i, j, xN_flag) - summer;
383 return line1 + line2 + line3 + line4;
384}
385CoolPropDbl MixtureDerivatives::d2_nd_ndalphardni_dnj_dDelta_dTau__constx(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
386 x_N_dependency_flag xN_flag) {
387 double line1 = d3_ndalphardni_dDelta2_dTau(HEOS, i, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag);
388 double line2 = d2_ndalphardni_dDelta_dTau(HEOS, i, xN_flag) * d_nddeltadni_dDelta(HEOS, j, xN_flag);
389 double line3 = d2_ndalphardni_dDelta_dTau(HEOS, i, xN_flag) * d_ndtaudni_dTau(HEOS, j, xN_flag);
390 double line4 = d3_ndalphardni_dDelta_dTau2(HEOS, i, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag);
391 double summer = 0;
392 std::size_t kmax = HEOS.mole_fractions.size();
393 if (xN_flag == XN_DEPENDENT) {
394 kmax--;
395 }
396 for (unsigned int k = 0; k < kmax; k++) {
397 summer += HEOS.mole_fractions[k] * d3_ndalphardni_dxj_dDelta_dTau__constxi(HEOS, i, k, xN_flag);
398 }
399 double line5 = d3_ndalphardni_dxj_dDelta_dTau__constxi(HEOS, i, j, xN_flag) - summer;
400 return line1 + line2 + line3 + line4 + line5;
401}
402CoolPropDbl MixtureDerivatives::d_nd_ndalphardni_dnj_dDelta__consttau_x(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
403 x_N_dependency_flag xN_flag) {
404 double line1 = d2_ndalphardni_dDelta2(HEOS, i, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag)
405 + d_ndalphardni_dDelta(HEOS, i, xN_flag) * d_nddeltadni_dDelta(HEOS, j, xN_flag);
406 double line2 = d2_ndalphardni_dDelta_dTau(HEOS, i, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag);
407 double summer = 0;
408 std::size_t kmax = HEOS.mole_fractions.size();
409 if (xN_flag == XN_DEPENDENT) {
410 kmax--;
411 }
412 for (unsigned int k = 0; k < kmax; k++) {
413 summer += HEOS.mole_fractions[k] * d2_ndalphardni_dxj_dDelta__consttau_xi(HEOS, i, k, xN_flag);
414 }
415 double line3 = d2_ndalphardni_dxj_dDelta__consttau_xi(HEOS, i, j, xN_flag) - summer;
416 return line1 + line2 + line3;
417}
418
419CoolPropDbl MixtureDerivatives::d_nd_ndalphardni_dnj_dxk__consttau_delta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
420 std::size_t k, x_N_dependency_flag xN_flag) {
421 double line1 = d_ndalphardni_dDelta(HEOS, i, xN_flag) * d_nddeltadni_dxj__constdelta_tau(HEOS, j, k, xN_flag)
422 + d2_ndalphardni_dxj_dDelta__consttau_xi(HEOS, i, k, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag);
423
424 double line2 = d_ndalphardni_dTau(HEOS, i, xN_flag) * d_ndtaudni_dxj__constdelta_tau(HEOS, j, k, xN_flag)
425 + d2_ndalphardni_dxj_dTau__constdelta_xi(HEOS, i, k, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag);
426
427 double line3 = d2_ndalphardni_dxj_dxk__constdelta_tau_xi(HEOS, i, j, k, xN_flag) - d_ndalphardni_dxj__constdelta_tau_xi(HEOS, i, k, xN_flag);
428
429 std::size_t mmax = HEOS.mole_fractions.size();
430 if (xN_flag == XN_DEPENDENT) {
431 mmax--;
432 }
433 for (unsigned int m = 0; m < mmax; m++) {
434 line3 -= HEOS.mole_fractions[m] * d2_ndalphardni_dxj_dxk__constdelta_tau_xi(HEOS, i, m, k, xN_flag);
435 }
436 return line1 + line2 + line3;
437}
438CoolPropDbl MixtureDerivatives::d2_nd_ndalphardni_dnj_dxk_dTau__constdelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
439 std::size_t k, x_N_dependency_flag xN_flag) {
440 double line1 = d2_ndalphardni_dDelta_dTau(HEOS, i, xN_flag) * d_nddeltadni_dxj__constdelta_tau(HEOS, j, k, xN_flag)
441 + d3_ndalphardni_dxj_dDelta_dTau__constxi(HEOS, i, k, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag);
442
443 double line2 = d_ndalphardni_dTau(HEOS, i, xN_flag) * d2_ndtaudni_dxj_dTau__constdelta(HEOS, j, k, xN_flag)
444 + d2_ndalphardni_dTau2(HEOS, i, xN_flag) * d_ndtaudni_dxj__constdelta_tau(HEOS, j, k, xN_flag)
445 + d2_ndalphardni_dxj_dTau__constdelta_xi(HEOS, i, k, xN_flag) * d_ndtaudni_dTau(HEOS, j, xN_flag)
446 + d3_ndalphardni_dxj_dTau2__constdelta_xi(HEOS, i, k, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag);
447
448 double line3 = d3_ndalphardni_dxj_dxk_dTau__constdelta_xi(HEOS, i, j, k, xN_flag) - d2_ndalphardni_dxj_dTau__constdelta_xi(HEOS, i, k, xN_flag);
449
450 std::size_t mmax = HEOS.mole_fractions.size();
451 if (xN_flag == XN_DEPENDENT) {
452 mmax--;
453 }
454 for (unsigned int m = 0; m < mmax; m++) {
455 line3 -= HEOS.mole_fractions[m] * d3_ndalphardni_dxj_dxk_dTau__constdelta_xi(HEOS, i, m, k, xN_flag);
456 }
457 return line1 + line2 + line3;
458}
459CoolPropDbl MixtureDerivatives::d2_nd_ndalphardni_dnj_dxk_dDelta__consttau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
460 std::size_t k, x_N_dependency_flag xN_flag) {
461 double line1 = d_ndalphardni_dDelta(HEOS, i, xN_flag) * d2_nddeltadni_dxj_dDelta__consttau(HEOS, j, k, xN_flag)
462 + d2_ndalphardni_dDelta2(HEOS, i, xN_flag) * d_nddeltadni_dxj__constdelta_tau(HEOS, j, k, xN_flag)
463 + d2_ndalphardni_dxj_dDelta__consttau_xi(HEOS, i, k, xN_flag) * d_nddeltadni_dDelta(HEOS, j, xN_flag)
464 + d3_ndalphardni_dxj_dDelta2__consttau_xi(HEOS, i, k, xN_flag) * nddeltadni__constT_V_nj(HEOS, j, xN_flag);
465
466 double line2 = d2_ndalphardni_dDelta_dTau(HEOS, i, xN_flag) * d_ndtaudni_dxj__constdelta_tau(HEOS, j, k, xN_flag)
467 + d3_ndalphardni_dxj_dDelta_dTau__constxi(HEOS, i, k, xN_flag) * ndtaudni__constT_V_nj(HEOS, j, xN_flag);
468
469 double line3 = d3_ndalphardni_dxj_dxk_dDelta__consttau_xi(HEOS, i, j, k, xN_flag) - d2_ndalphardni_dxj_dDelta__consttau_xi(HEOS, i, k, xN_flag);
470 std::size_t mmax = HEOS.mole_fractions.size();
471 if (xN_flag == XN_DEPENDENT) {
472 mmax--;
473 }
474 for (unsigned int m = 0; m < mmax; m++) {
475 line3 -= HEOS.mole_fractions[m] * d3_ndalphardni_dxj_dxk_dDelta__consttau_xi(HEOS, i, m, k, xN_flag);
476 }
477 return line1 + line2 + line3;
478}
479CoolPropDbl MixtureDerivatives::d_ndalphardni_dDelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
480 // The first line
481 double term1 = (HEOS._delta.pt() * HEOS.d2alphar_dDelta2() + HEOS.dalphar_dDelta())
482 * (1 - 1 / HEOS._reducing.rhomolar * HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag));
483
484 // The second line
485 double term2 =
486 HEOS._tau.pt() * HEOS.d2alphar_dDelta_dTau() * (1 / HEOS._reducing.T) * HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag);
487
488 // The third line
489 double term3 = HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, i, xN_flag);
490 std::size_t kmax = HEOS.mole_fractions.size();
491 if (xN_flag == XN_DEPENDENT) {
492 kmax--;
493 }
494 for (unsigned int k = 0; k < kmax; k++) {
495 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, k, xN_flag);
496 }
497 return term1 + term2 + term3;
498}
499CoolPropDbl MixtureDerivatives::d2_ndalphardni_dDelta2(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
500 double term1 = (2 * HEOS.d2alphar_dDelta2() + HEOS.delta() * HEOS.d3alphar_dDelta3()) * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
501 double term2 = HEOS.tau() * HEOS.d3alphar_dDelta2_dTau() * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
502 double term3 = HEOS.residual_helmholtz->d3alphar_dxi_dDelta2(HEOS, i, xN_flag);
503 std::size_t kmax = HEOS.mole_fractions.size();
504 if (xN_flag == XN_DEPENDENT) {
505 kmax--;
506 }
507 for (unsigned int k = 0; k < kmax; k++) {
508 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d3alphar_dxi_dDelta2(HEOS, k, xN_flag);
509 }
510 return term1 + term2 + term3;
511}
512CoolPropDbl MixtureDerivatives::d3_ndalphardni_dDelta3(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
513 double term1 = (3 * HEOS.d3alphar_dDelta3() + HEOS.delta() * HEOS.d4alphar_dDelta4()) * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
514 double term2 = HEOS.tau() * HEOS.d4alphar_dDelta3_dTau() * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
515 double term3 = HEOS.residual_helmholtz->d4alphar_dxi_dDelta3(HEOS, i, xN_flag);
516 std::size_t kmax = HEOS.mole_fractions.size();
517 if (xN_flag == XN_DEPENDENT) {
518 kmax--;
519 }
520 for (unsigned int k = 0; k < kmax; k++) {
521 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d4alphar_dxi_dDelta3(HEOS, k, xN_flag);
522 }
523 return term1 + term2 + term3;
524}
525CoolPropDbl MixtureDerivatives::d2_ndalphardni_dDelta_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
526 double term1 =
527 (HEOS.d2alphar_dDelta_dTau() + HEOS.delta() * HEOS.d3alphar_dDelta2_dTau()) * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
528 double term2 = (HEOS.tau() * HEOS.d3alphar_dDelta_dTau2() + HEOS.d2alphar_dDelta_dTau()) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
529 double term3 = HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, i, xN_flag);
530 std::size_t kmax = HEOS.mole_fractions.size();
531 if (xN_flag == XN_DEPENDENT) {
532 kmax--;
533 }
534 for (unsigned int k = 0; k < kmax; k++) {
535 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, k, xN_flag);
536 }
537 return term1 + term2 + term3;
538}
539CoolPropDbl MixtureDerivatives::d3_ndalphardni_dDelta2_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
540 double term1 =
541 (2 * HEOS.d3alphar_dDelta2_dTau() + HEOS.delta() * HEOS.d4alphar_dDelta3_dTau()) * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
542 double term2 =
543 (HEOS.tau() * HEOS.d4alphar_dDelta2_dTau2() + HEOS.d3alphar_dDelta2_dTau()) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
544 double term3 = HEOS.residual_helmholtz->d4alphar_dxi_dDelta2_dTau(HEOS, i, xN_flag);
545 std::size_t kmax = HEOS.mole_fractions.size();
546 if (xN_flag == XN_DEPENDENT) {
547 kmax--;
548 }
549 for (unsigned int k = 0; k < kmax; k++) {
550 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d4alphar_dxi_dDelta2_dTau(HEOS, k, xN_flag);
551 }
552 return term1 + term2 + term3;
553}
554CoolPropDbl MixtureDerivatives::d3_ndalphardni_dDelta_dTau2(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
555 double term1 =
556 (HEOS.d3alphar_dDelta_dTau2() + HEOS.delta() * HEOS.d4alphar_dDelta2_dTau2()) * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
557 double term2 =
558 (HEOS.tau() * HEOS.d4alphar_dDelta_dTau3() + 2 * HEOS.d3alphar_dDelta_dTau2()) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
559 double term3 = HEOS.residual_helmholtz->d4alphar_dxi_dDelta_dTau2(HEOS, i, xN_flag);
560 std::size_t kmax = HEOS.mole_fractions.size();
561 if (xN_flag == XN_DEPENDENT) {
562 kmax--;
563 }
564 for (unsigned int k = 0; k < kmax; k++) {
565 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d4alphar_dxi_dDelta_dTau2(HEOS, k, xN_flag);
566 }
567 return term1 + term2 + term3;
568}
569CoolPropDbl MixtureDerivatives::d2_ndalphardni_dTau2(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
570 double term1 = HEOS.delta() * HEOS.d3alphar_dDelta_dTau2() * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
571 double term2 = (2 * HEOS.d2alphar_dTau2() + HEOS.tau() * HEOS.d3alphar_dTau3()) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
572 double term3 = HEOS.residual_helmholtz->d3alphar_dxi_dTau2(HEOS, i, xN_flag);
573 std::size_t kmax = HEOS.mole_fractions.size();
574 if (xN_flag == XN_DEPENDENT) {
575 kmax--;
576 }
577 for (unsigned int k = 0; k < kmax; k++) {
578 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d3alphar_dxi_dTau2(HEOS, k, xN_flag);
579 }
580 return term1 + term2 + term3;
581}
582CoolPropDbl MixtureDerivatives::d3_ndalphardni_dTau3(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
583 double term1 = HEOS.delta() * HEOS.d4alphar_dDelta_dTau3() * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
584 double term2 = (3 * HEOS.d3alphar_dTau3() + HEOS.tau() * HEOS.d4alphar_dTau4()) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
585 double term3 = HEOS.residual_helmholtz->d4alphar_dxi_dTau3(HEOS, i, xN_flag);
586 std::size_t kmax = HEOS.mole_fractions.size();
587 if (xN_flag == XN_DEPENDENT) {
588 kmax--;
589 }
590 for (unsigned int k = 0; k < kmax; k++) {
591 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d4alphar_dxi_dTau3(HEOS, k, xN_flag);
592 }
593 return term1 + term2 + term3;
594}
595
596CoolPropDbl MixtureDerivatives::d2_ndalphardni_dxj_dDelta__consttau_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
597 x_N_dependency_flag xN_flag) {
598 double term1 =
599 (HEOS.dalphar_dDelta() + HEOS.delta() * HEOS.d2alphar_dDelta2()) * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, j, xN_flag);
600 double term2 = (HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, j, xN_flag)
601 + HEOS.delta() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta2(HEOS, j, xN_flag))
602 * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
603 double term3 = HEOS.tau() * HEOS.d2alphar_dDelta_dTau() * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, j, xN_flag);
604 double term4 =
605 HEOS.tau() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, j, xN_flag) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
606 double term5 = HEOS.residual_helmholtz->d3alphar_dxi_dxj_dDelta(HEOS, i, j, xN_flag);
607 std::size_t kmax = HEOS.mole_fractions.size();
608 if (xN_flag == XN_DEPENDENT) {
609 kmax--;
610 }
611 for (unsigned int k = 0; k < kmax; k++) {
612 term5 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d3alphar_dxi_dxj_dDelta(HEOS, k, j, xN_flag)
613 + Kronecker_delta(k, j) * HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, k, xN_flag);
614 }
615 return term1 + term2 + term3 + term4 + term5;
616}
617CoolPropDbl MixtureDerivatives::d2_ndalphardni_dxj_dTau__constdelta_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
618 x_N_dependency_flag xN_flag) {
619 double term1 = HEOS.delta() * HEOS.d2alphar_dDelta_dTau() * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, j, xN_flag);
620 double term2 =
621 HEOS.delta() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, j, xN_flag) * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
622 double term3 = (HEOS.tau() * HEOS.d2alphar_dTau2() + HEOS.dalphar_dTau()) * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, j, xN_flag);
623 double term4 =
624 (HEOS.tau() * HEOS.residual_helmholtz->d3alphar_dxi_dTau2(HEOS, j, xN_flag) + HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, j, xN_flag))
625 * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
626 double term5 = HEOS.residual_helmholtz->d3alphar_dxi_dxj_dTau(HEOS, i, j, xN_flag);
627 std::size_t kmax = HEOS.mole_fractions.size();
628 if (xN_flag == XN_DEPENDENT) {
629 kmax--;
630 }
631 for (unsigned int k = 0; k < kmax; k++) {
632 term5 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d3alphar_dxi_dxj_dTau(HEOS, k, j, xN_flag)
633 + Kronecker_delta(k, j) * HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, k, xN_flag);
634 }
635 return term1 + term2 + term3 + term4 + term5;
636}
637CoolPropDbl MixtureDerivatives::d3_ndalphardni_dxj_dTau2__constdelta_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
638 x_N_dependency_flag xN_flag) {
639 double term1 = HEOS.delta() * HEOS.d3alphar_dDelta_dTau2() * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, j, xN_flag);
640 double term2 =
641 HEOS.delta() * HEOS.residual_helmholtz->d4alphar_dxi_dDelta_dTau2(HEOS, j, xN_flag) * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
642 double term3 = (HEOS.tau() * HEOS.d3alphar_dTau3() + 2 * HEOS.d2alphar_dTau2()) * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, j, xN_flag);
643 double term4 =
644 (HEOS.tau() * HEOS.residual_helmholtz->d4alphar_dxi_dTau3(HEOS, j, xN_flag) + 2 * HEOS.residual_helmholtz->d3alphar_dxi_dTau2(HEOS, j, xN_flag))
645 * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
646 double term5 = HEOS.residual_helmholtz->d4alphar_dxi_dxj_dTau2(HEOS, i, j, xN_flag);
647 std::size_t kmax = HEOS.mole_fractions.size();
648 if (xN_flag == XN_DEPENDENT) {
649 kmax--;
650 }
651 for (unsigned int k = 0; k < kmax; k++) {
652 term5 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dTau2(HEOS, k, j, xN_flag)
653 + Kronecker_delta(k, j) * HEOS.residual_helmholtz->d3alphar_dxi_dTau2(HEOS, k, xN_flag);
654 }
655 return term1 + term2 + term3 + term4 + term5;
656}
657CoolPropDbl MixtureDerivatives::d3_ndalphardni_dxj_dDelta2__consttau_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
658 x_N_dependency_flag xN_flag) {
659 double term1 =
660 (2 * HEOS.d2alphar_dDelta2() + HEOS.delta() * HEOS.d3alphar_dDelta3()) * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, j, xN_flag);
661 double term2 = (2 * HEOS.residual_helmholtz->d3alphar_dxi_dDelta2(HEOS, j, xN_flag)
662 + HEOS.delta() * HEOS.residual_helmholtz->d4alphar_dxi_dDelta3(HEOS, j, xN_flag))
663 * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
664 double term3 = HEOS.tau() * HEOS.d3alphar_dDelta2_dTau() * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, j, xN_flag);
665 double term4 =
666 HEOS.tau() * HEOS.residual_helmholtz->d4alphar_dxi_dDelta2_dTau(HEOS, j, xN_flag) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
667 double term5 = HEOS.residual_helmholtz->d4alphar_dxi_dxj_dDelta2(HEOS, i, j, xN_flag);
668 std::size_t kmax = HEOS.mole_fractions.size();
669 if (xN_flag == XN_DEPENDENT) {
670 kmax--;
671 }
672 for (unsigned int k = 0; k < kmax; k++) {
673 term5 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dDelta2(HEOS, k, j, xN_flag)
674 + Kronecker_delta(k, j) * HEOS.residual_helmholtz->d3alphar_dxi_dDelta2(HEOS, k, xN_flag);
675 }
676 return term1 + term2 + term3 + term4 + term5;
677}
678CoolPropDbl MixtureDerivatives::d3_ndalphardni_dxj_dDelta_dTau__constxi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
679 x_N_dependency_flag xN_flag) {
680 double term1 =
681 (HEOS.d2alphar_dDelta_dTau() + HEOS.delta() * HEOS.d3alphar_dDelta2_dTau()) * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, j, xN_flag);
682 double term2 = (HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, j, xN_flag)
683 + HEOS.delta() * HEOS.residual_helmholtz->d4alphar_dxi_dDelta2_dTau(HEOS, j, xN_flag))
684 * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
685 double term3 =
686 (HEOS.tau() * HEOS.d3alphar_dDelta_dTau2() + HEOS.d2alphar_dDelta_dTau()) * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, j, xN_flag);
687 double term4 = (HEOS.tau() * HEOS.residual_helmholtz->d4alphar_dxi_dDelta_dTau2(HEOS, j, xN_flag)
688 + HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, j, xN_flag))
689 * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
690 double term5 = HEOS.residual_helmholtz->d4alphar_dxi_dxj_dDelta_dTau(HEOS, i, j, xN_flag);
691 std::size_t kmax = HEOS.mole_fractions.size();
692 if (xN_flag == XN_DEPENDENT) {
693 kmax--;
694 }
695 for (unsigned int k = 0; k < kmax; k++) {
696 term5 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dDelta_dTau(HEOS, k, j, xN_flag)
697 + Kronecker_delta(k, j) * HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, k, xN_flag);
698 }
699 return term1 + term2 + term3 + term4 + term5;
700}
701
702CoolPropDbl MixtureDerivatives::d_ndalphardni_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
703 // The first line
704 double term1 = HEOS._delta.pt() * HEOS.d2alphar_dDelta_dTau()
705 * (1 - 1 / HEOS._reducing.rhomolar * HEOS.Reducing->ndrhorbardni__constnj(HEOS.mole_fractions, i, xN_flag));
706
707 // The second line
708 double term2 = (HEOS._tau.pt() * HEOS.d2alphar_dTau2() + HEOS.dalphar_dTau()) * (1 / HEOS._reducing.T)
709 * HEOS.Reducing->ndTrdni__constnj(HEOS.mole_fractions, i, xN_flag);
710
711 // The third line
712 double term3 = HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, i, xN_flag);
713 std::size_t kmax = HEOS.mole_fractions.size();
714 if (xN_flag == XN_DEPENDENT) {
715 kmax--;
716 }
717 for (unsigned int k = 0; k < kmax; k++) {
718 term3 -= HEOS.mole_fractions[k] * HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, k, xN_flag);
719 }
720 return term1 + term2 + term3;
721}
722
723CoolPropDbl MixtureDerivatives::d2_ndalphardni_dxj_dxk__constdelta_tau_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
724 std::size_t k, x_N_dependency_flag xN_flag) {
725 double term1 =
726 HEOS.delta()
727 * (HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, j, xN_flag) * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, k, xN_flag)
728 + HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, k, xN_flag) * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, j, xN_flag));
729 double term2 =
730 HEOS.delta() * HEOS.residual_helmholtz->d3alphar_dxi_dxj_dDelta(HEOS, j, k, xN_flag) * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
731 double term3 = HEOS.delta() * HEOS.dalphar_dDelta() * HEOS.Reducing->d2_PSI_rho_dxj_dxk(HEOS.mole_fractions, i, j, k, xN_flag);
732
733 double term4 =
734 HEOS.tau()
735 * (HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, j, xN_flag) * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, k, xN_flag)
736 + HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, k, xN_flag) * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, j, xN_flag));
737 double term5 =
738 HEOS.tau() * HEOS.residual_helmholtz->d3alphar_dxi_dxj_dTau(HEOS, j, k, xN_flag) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
739 double term6 = HEOS.tau() * HEOS.dalphar_dTau() * HEOS.Reducing->d2_PSI_T_dxj_dxk(HEOS.mole_fractions, i, j, k, xN_flag);
740
741 double term7 =
742 HEOS.residual_helmholtz->d3alphardxidxjdxk(HEOS, i, j, k, xN_flag) - 2 * HEOS.residual_helmholtz->d2alphardxidxj(HEOS, j, k, xN_flag);
743 std::size_t mmax = HEOS.mole_fractions.size();
744 if (xN_flag == XN_DEPENDENT) {
745 mmax--;
746 }
747 for (unsigned int m = 0; m < mmax; m++) {
748 term7 -= HEOS.mole_fractions[m] * HEOS.residual_helmholtz->d3alphardxidxjdxk(HEOS, j, k, m, xN_flag);
749 }
750
751 return term1 + term2 + term3 + term4 + term5 + term6 + term7;
752}
753
754CoolPropDbl MixtureDerivatives::d3_ndalphardni_dxj_dxk_dTau__constdelta_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
755 std::size_t k, x_N_dependency_flag xN_flag) {
756 double term1a = HEOS.delta() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, j, xN_flag)
757 * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, k, xN_flag);
758 double term1b = HEOS.delta() * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dDelta_dTau(HEOS, j, k, xN_flag)
759 * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
760 double term1c = HEOS.delta() * HEOS.d2alphar_dDelta_dTau() * HEOS.Reducing->d2_PSI_rho_dxj_dxk(HEOS.mole_fractions, i, j, k, xN_flag);
761 double term1d = HEOS.delta() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, k, xN_flag)
762 * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, j, xN_flag);
763 double term1 = term1a + term1b + term1c + term1d;
764
765 double term2a =
766 (HEOS.tau() * HEOS.residual_helmholtz->d3alphar_dxi_dTau2(HEOS, j, xN_flag) + HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, j, xN_flag))
767 * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, k, xN_flag);
768 double term2b = (HEOS.tau() * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dTau2(HEOS, j, k, xN_flag)
769 + HEOS.residual_helmholtz->d3alphar_dxi_dxj_dTau(HEOS, j, k, xN_flag))
770 * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
771 double term2c =
772 (HEOS.tau() * HEOS.d2alphar_dTau2() + HEOS.dalphar_dTau()) * HEOS.Reducing->d2_PSI_T_dxj_dxk(HEOS.mole_fractions, i, j, k, xN_flag);
773 double term2d =
774 (HEOS.tau() * HEOS.residual_helmholtz->d3alphar_dxi_dTau2(HEOS, k, xN_flag) + HEOS.residual_helmholtz->d2alphar_dxi_dTau(HEOS, k, xN_flag))
775 * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, j, xN_flag);
776 double term2 = term2a + term2b + term2c + term2d;
777
778 double term3 = HEOS.residual_helmholtz->d4alphar_dxi_dxj_dxk_dTau(HEOS, i, j, k, xN_flag)
779 - 2 * HEOS.residual_helmholtz->d3alphar_dxi_dxj_dTau(HEOS, j, k, xN_flag);
780 std::size_t mmax = HEOS.mole_fractions.size();
781 if (xN_flag == XN_DEPENDENT) {
782 mmax--;
783 }
784 for (unsigned int m = 0; m < mmax; m++) {
785 term3 -= HEOS.mole_fractions[m] * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dxk_dTau(HEOS, j, k, m, xN_flag);
786 }
787
788 return term1 + term2 + term3;
789}
790
791CoolPropDbl MixtureDerivatives::d3_ndalphardni_dxj_dxk_dDelta__consttau_xi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
792 std::size_t k, x_N_dependency_flag xN_flag) {
793 double term1a = (HEOS.delta() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta2(HEOS, j, xN_flag)
794 + HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, j, xN_flag))
795 * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, k, xN_flag);
796 double term1b = (HEOS.delta() * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dDelta2(HEOS, j, k, xN_flag)
797 + HEOS.residual_helmholtz->d3alphar_dxi_dxj_dDelta(HEOS, j, k, xN_flag))
798 * HEOS.Reducing->PSI_rho(HEOS.mole_fractions, i, xN_flag);
799 double term1c =
800 (HEOS.delta() * HEOS.d2alphar_dDelta2() + HEOS.dalphar_dDelta()) * HEOS.Reducing->d2_PSI_rho_dxj_dxk(HEOS.mole_fractions, i, j, k, xN_flag);
801 double term1d = (HEOS.delta() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta2(HEOS, k, xN_flag)
802 + HEOS.residual_helmholtz->d2alphar_dxi_dDelta(HEOS, k, xN_flag))
803 * HEOS.Reducing->d_PSI_rho_dxj(HEOS.mole_fractions, i, j, xN_flag);
804 double term1 = term1a + term1b + term1c + term1d;
805
806 double term2a = HEOS.tau() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, j, xN_flag)
807 * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, k, xN_flag);
808 double term2b =
809 HEOS.tau() * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dDelta_dTau(HEOS, j, k, xN_flag) * HEOS.Reducing->PSI_T(HEOS.mole_fractions, i, xN_flag);
810 double term2c = HEOS.tau() * HEOS.d2alphar_dDelta_dTau() * HEOS.Reducing->d2_PSI_T_dxj_dxk(HEOS.mole_fractions, i, j, k, xN_flag);
811 double term2d = HEOS.tau() * HEOS.residual_helmholtz->d3alphar_dxi_dDelta_dTau(HEOS, k, xN_flag)
812 * HEOS.Reducing->d_PSI_T_dxj(HEOS.mole_fractions, i, j, xN_flag);
813 double term2 = term2a + term2b + term2c + term2d;
814
815 double term3 = HEOS.residual_helmholtz->d4alphar_dxi_dxj_dxk_dDelta(HEOS, i, j, k, xN_flag)
816 - 2 * HEOS.residual_helmholtz->d3alphar_dxi_dxj_dDelta(HEOS, j, k, xN_flag);
817 std::size_t mmax = HEOS.mole_fractions.size();
818 if (xN_flag == XN_DEPENDENT) {
819 mmax--;
820 }
821 for (unsigned int m = 0; m < mmax; m++) {
822 term3 -= HEOS.mole_fractions[m] * HEOS.residual_helmholtz->d4alphar_dxi_dxj_dxk_dDelta(HEOS, j, k, m, xN_flag);
823 }
824 return term1 + term2 + term3;
825}
826
827CoolPropDbl MixtureDerivatives::dalpha0_dxi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
828 // Reducing values are constant for all components under consideration
829 double Tr = HEOS.T_reducing();
830 double rhor = HEOS.rhomolar_reducing();
831
832 // Values for the i-th component
833 double Tci = HEOS.get_fluid_constant(i, iT_critical);
834 double rhoci = HEOS.get_fluid_constant(i, irhomolar_critical);
835 double tau_oi = HEOS.tau() * Tci / Tr;
836 double delta_oi = HEOS.delta() * rhor / rhoci;
837 double Rratioi = 1; //HEOS.gas_constant()/HEOS.components[i].EOS().R_u;
838
839 double logxi = (std::abs(HEOS.mole_fractions[i]) > DBL_EPSILON) ? (log(HEOS.mole_fractions[i])) : 0;
840 double term = Rratioi * HEOS.components[i].EOS().alpha0.base(tau_oi, delta_oi) + logxi + 1;
841
842 std::size_t kmax = HEOS.mole_fractions.size();
843 if (xN_flag == XN_DEPENDENT) {
844 kmax--;
845 }
846 for (std::size_t k = 0; k < kmax; ++k) {
847 double xk = HEOS.mole_fractions[k];
848 double Tck = HEOS.get_fluid_constant(k, iT_critical);
849 double rhock = HEOS.get_fluid_constant(k, irhomolar_critical);
850 double tau_ok = HEOS.tau() * Tck / Tr;
851 double delta_ok = HEOS.delta() * rhor / rhock;
852 double dtauok_dxi = -tau_ok / Tr * HEOS.Reducing->dTrdxi__constxj(HEOS.mole_fractions, i, xN_flag); // (Gernert, supp, B.19)
853 double ddeltaok_dxi = delta_ok / rhor * HEOS.Reducing->drhormolardxi__constxj(HEOS.mole_fractions, i, xN_flag); // (Gernert, supp. B.20)
854
855 double Rratiok = 1; //HEOS.gas_constant()/HEOS.components[k].EOS().R_u;
856 HelmholtzDerivatives alpha0kterms = HEOS.components[k].EOS().alpha0.all(tau_ok, delta_ok);
857 double dalpha0_ok_dxi = alpha0kterms.dalphar_dtau * dtauok_dxi + alpha0kterms.dalphar_ddelta * ddeltaok_dxi;
858 term += xk * (Rratiok * dalpha0_ok_dxi);
859 }
860 return term;
861}
862
863CoolPropDbl MixtureDerivatives::d2alpha0_dxi_dDelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
864 // Reducing values are constant for all components under consideration
865 double Tr = HEOS.T_reducing();
866 double rhor = HEOS.rhomolar_reducing();
867
868 // Values for the i-th component
869 double Tci = HEOS.get_fluid_constant(i, iT_critical);
870 double rhoci = HEOS.get_fluid_constant(i, irhomolar_critical);
871 double tau_oi = HEOS.tau() * Tci / Tr;
872 double delta_oi = HEOS.delta() * rhor / rhoci;
873 double Rratioi = 1; //HEOS.gas_constant()/HEOS.components[i].EOS().R_u;
874
875 double term = rhor / rhoci * Rratioi * HEOS.components[i].EOS().alpha0.dDelta(tau_oi, delta_oi);
876
877 std::size_t kmax = HEOS.mole_fractions.size();
878 if (xN_flag == XN_DEPENDENT) {
879 kmax--;
880 }
881 for (std::size_t k = 0; k < kmax; ++k) {
882 double xk = HEOS.mole_fractions[k];
883 double Tck = HEOS.get_fluid_constant(k, iT_critical);
884 double rhock = HEOS.get_fluid_constant(k, irhomolar_critical);
885 double tau_ok = HEOS.tau() * Tck / Tr;
886 double delta_ok = HEOS.delta() * rhor / rhock;
887 double dtauok_dxi = -tau_ok / Tr * HEOS.Reducing->dTrdxi__constxj(HEOS.mole_fractions, i, xN_flag); // (Gernert, supp, B.19)
888 double drhor_dxi = HEOS.Reducing->drhormolardxi__constxj(HEOS.mole_fractions, i, xN_flag);
889 double ddeltaok_dxi = delta_ok / rhor * drhor_dxi; // (Gernert, supp. B.20)
890
891 //double Rratiok = 1;//HEOS.gas_constant()/HEOS.components[k].EOS().R_u;
892 HelmholtzDerivatives alpha0kterms = HEOS.components[k].EOS().alpha0.all(tau_ok, delta_ok);
893 double dalpha0ok_ddeltaok = alpha0kterms.dalphar_ddelta;
894
895 double d_dalpha0ok_ddeltaok_dxi = alpha0kterms.d2alphar_ddelta_dtau * dtauok_dxi + alpha0kterms.d2alphar_ddelta2 * ddeltaok_dxi;
896 term += xk / rhock * (rhor * d_dalpha0ok_ddeltaok_dxi + drhor_dxi * dalpha0ok_ddeltaok);
897 }
898 return term;
899}
900
901CoolPropDbl MixtureDerivatives::d2alpha0_dxi_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
902 // Reducing values are constant for all components under consideration
903 double Tr = HEOS.T_reducing();
904 double rhor = HEOS.rhomolar_reducing();
905
906 // Values for the i-th component
907 double Tci = HEOS.get_fluid_constant(i, iT_critical);
908 double rhoci = HEOS.get_fluid_constant(i, irhomolar_critical);
909 double tau_oi = HEOS.tau() * Tci / Tr;
910 double delta_oi = HEOS.delta() * rhor / rhoci;
911 double Rratioi = 1; //HEOS.gas_constant()/HEOS.components[i].EOS().R_u;
912
913 double term = Tci / Tr * Rratioi * HEOS.components[i].EOS().alpha0.dTau(tau_oi, delta_oi);
914
915 std::size_t kmax = HEOS.mole_fractions.size();
916 if (xN_flag == XN_DEPENDENT) {
917 kmax--;
918 }
919 for (std::size_t k = 0; k < kmax; ++k) {
920 double xk = HEOS.mole_fractions[k];
921 double Tck = HEOS.get_fluid_constant(k, iT_critical);
922 double rhock = HEOS.get_fluid_constant(k, irhomolar_critical);
923 double tau_ok = HEOS.tau() * Tck / Tr;
924 double delta_ok = HEOS.delta() * rhor / rhock;
925 double dTr_dxi = HEOS.Reducing->dTrdxi__constxj(HEOS.mole_fractions, i, xN_flag);
926 double dtauok_dxi = -tau_ok / Tr * dTr_dxi; // (Gernert, supp, B.19)
927 double drhor_dxi = HEOS.Reducing->drhormolardxi__constxj(HEOS.mole_fractions, i, xN_flag);
928 double ddeltaok_dxi = delta_ok / rhor * drhor_dxi; // (Gernert, supp. B.20)
929
930 //double Rratiok = 1;//HEOS.gas_constant()/HEOS.components[k].EOS().R_u;
931 HelmholtzDerivatives alpha0kterms = HEOS.components[k].EOS().alpha0.all(tau_ok, delta_ok);
932 double dalpha0ok_dtauok = alpha0kterms.dalphar_dtau;
933 double d_dalpha0ok_dTauok_dxi = alpha0kterms.d2alphar_dtau2 * dtauok_dxi + alpha0kterms.d2alphar_ddelta_dtau * ddeltaok_dxi;
934 term += xk * Tck * (1 / Tr * d_dalpha0ok_dTauok_dxi + -1 / POW2(Tr) * dTr_dxi * dalpha0ok_dtauok);
935 }
936 return term;
937}
938
939CoolPropDbl MixtureDerivatives::d2alpha0dxidxj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
940 // Reducing values are constant for all components under consideration
941 double Tr = HEOS.T_reducing();
942 double rhor = HEOS.rhomolar_reducing();
943
944 // Values for the i-th component
945 double Tci = HEOS.get_fluid_constant(i, iT_critical);
946 double rhoci = HEOS.get_fluid_constant(i, irhomolar_critical);
947 double tau_oi = HEOS.tau() * Tci / Tr;
948 double delta_oi = HEOS.delta() * rhor / rhoci;
949 double dTr_dxi = HEOS.Reducing->dTrdxi__constxj(HEOS.mole_fractions, i, xN_flag);
950 double drhor_dxi = HEOS.Reducing->drhormolardxi__constxj(HEOS.mole_fractions, i, xN_flag);
951
952 // Values for the j-th component
953 double Tcj = HEOS.get_fluid_constant(j, iT_critical);
954 double rhocj = HEOS.get_fluid_constant(j, irhomolar_critical);
955 double tau_oj = HEOS.tau() * Tcj / Tr;
956 double delta_oj = HEOS.delta() * rhor / rhocj;
957 double dTr_dxj = HEOS.Reducing->dTrdxi__constxj(HEOS.mole_fractions, j, xN_flag);
958 double drhor_dxj = HEOS.Reducing->drhormolardxi__constxj(HEOS.mole_fractions, j, xN_flag);
959
960 // Cross-terms with both i & j
961 double dtauoi_dxj = -tau_oi / Tr * dTr_dxj; // (Gernert, supp, B.19)
962 double ddeltaoi_dxj = delta_oi / rhor * drhor_dxj; // (Gernert, supp. B.20)
963 double dtauoj_dxi = -tau_oj / Tr * dTr_dxi; // (Gernert, supp, B.19)
964 double ddeltaoj_dxi = delta_oj / rhor * drhor_dxi; // (Gernert, supp. B.20)
965 double d2Tr_dxidxj = HEOS.Reducing->d2Trdxidxj(HEOS.mole_fractions, i, j, xN_flag);
966 double d2rhor_dxidxj = HEOS.Reducing->d2rhormolardxidxj(HEOS.mole_fractions, i, j, xN_flag);
967
968 //double Rratioi = 1;//HEOS.gas_constant()/HEOS.components[i].EOS().R_u;
969 HelmholtzDerivatives alpha0iterms = HEOS.components[i].EOS().alpha0.all(tau_oi, delta_oi),
970 alpha0jterms = HEOS.components[j].EOS().alpha0.all(tau_oj, delta_oj);
971
972 double d_dalpha0oi_dxj = alpha0iterms.dalphar_dtau * dtauoi_dxj + alpha0iterms.dalphar_ddelta * ddeltaoi_dxj;
973 double d_dalpha0oj_dxi = alpha0jterms.dalphar_dtau * dtauoj_dxi + alpha0jterms.dalphar_ddelta * ddeltaoj_dxi;
974
975 double xi = HEOS.mole_fractions[i], xj = HEOS.mole_fractions[j];
976 double Kronecker_delta_over_xi = (xj > DBL_EPSILON && xi > DBL_EPSILON) ? Kronecker_delta(i, j) / xi : 0;
977 double term = d_dalpha0oi_dxj + d_dalpha0oj_dxi + Kronecker_delta_over_xi;
978
979 std::size_t kmax = HEOS.mole_fractions.size();
980 if (xN_flag == XN_DEPENDENT) {
981 kmax--;
982 }
983 for (std::size_t k = 0; k < kmax; ++k) {
984 // Values for the k-th component
985 double xk = HEOS.mole_fractions[k];
986 double Tck = HEOS.get_fluid_constant(k, iT_critical);
987 double rhock = HEOS.get_fluid_constant(k, irhomolar_critical);
988 double tau_ok = HEOS.tau() * Tck / Tr;
989 double delta_ok = HEOS.delta() * rhor / rhock;
990
991 double dtauok_dxj = -tau_ok / Tr * dTr_dxj; // (Gernert, supp, B.19)
992 double ddeltaok_dxj = delta_ok / rhor * drhor_dxj; // (Gernert, supp. B.20)
993 double dtauok_dxi = -tau_ok / Tr * dTr_dxi; // (Gernert, supp, B.19)
994 double ddeltaok_dxi = delta_ok / rhor * drhor_dxi; // (Gernert, supp. B.20)
995
996 HelmholtzDerivatives alpha0kterms = HEOS.components[k].EOS().alpha0.all(tau_ok, delta_ok);
997 double dalpha0ok_dtauok = alpha0kterms.dalphar_dtau;
998 double d2tauok_dxidxj = -Tck * HEOS.tau() * (POW2(Tr) * d2Tr_dxidxj - dTr_dxi * (2 * Tr * dTr_dxj)) / POW4(Tr);
999 double d_dalpha0ok_dtauok_dxj = alpha0kterms.d2alphar_dtau2 * dtauok_dxj + alpha0kterms.d2alphar_ddelta_dtau * ddeltaok_dxj;
1000
1001 double dalpha0ok_ddeltaok = alpha0kterms.dalphar_ddelta;
1002 double d2deltaok_dxidxj = HEOS.delta() / rhock * d2rhor_dxidxj;
1003 double d_dalpha0ok_ddeltaok_dxj = alpha0kterms.d2alphar_ddelta_dtau * dtauok_dxj + alpha0kterms.d2alphar_ddelta2 * ddeltaok_dxj;
1004
1005 term += xk
1006 * (dalpha0ok_dtauok * d2tauok_dxidxj + d_dalpha0ok_dtauok_dxj * dtauok_dxi + dalpha0ok_ddeltaok * d2deltaok_dxidxj
1007 + d_dalpha0ok_ddeltaok_dxj * ddeltaok_dxi);
1008 }
1009 return term;
1010}
1011
1013CoolPropDbl MixtureDerivatives::dpsi_dDelta(HelmholtzEOSMixtureBackend& HEOS, x_N_dependency_flag xN_flag) {
1014 return HEOS.rhomolar_reducing() * HEOS.gas_constant() * HEOS.T() * (HEOS.delta() * dalpha_dDelta(HEOS) + alpha(HEOS));
1015}
1016
1018CoolPropDbl MixtureDerivatives::dpsi_dTau(HelmholtzEOSMixtureBackend& HEOS, x_N_dependency_flag xN_flag) {
1019 return HEOS.rhomolar_reducing() * HEOS.delta() * HEOS.gas_constant() * HEOS.T() / HEOS.tau() * (HEOS.tau() * dalpha_dTau(HEOS) - alpha(HEOS));
1020}
1021
1023CoolPropDbl MixtureDerivatives::dpsir_dTau(HelmholtzEOSMixtureBackend& HEOS, x_N_dependency_flag xN_flag) {
1024 return HEOS.rhomolar_reducing() * HEOS.delta() * HEOS.gas_constant() * HEOS.T() / HEOS.tau() * (HEOS.tau() * dalphar_dTau(HEOS) - alphar(HEOS));
1025}
1026
1027CoolPropDbl MixtureDerivatives::dpsi_dxi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
1028 return HEOS.delta() * HEOS.gas_constant() / HEOS.tau()
1029 * (alpha(HEOS, xN_flag) * d_rhorTr_dxi(HEOS, i, xN_flag) + HEOS.rhomolar_reducing() * HEOS.T_reducing() * dalpha_dxi(HEOS, i, xN_flag));
1030}
1031
1032CoolPropDbl MixtureDerivatives::dpsir_dxi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
1033 return HEOS.delta() * HEOS.gas_constant() / HEOS.tau()
1034 * (alphar(HEOS, xN_flag) * d_rhorTr_dxi(HEOS, i, xN_flag) + HEOS.rhomolar_reducing() * HEOS.T_reducing() * dalphar_dxi(HEOS, i, xN_flag));
1035}
1036
1037CoolPropDbl MixtureDerivatives::d_rhorTr_dxi(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
1038 GERG2008ReducingFunction* GERG = static_cast<GERG2008ReducingFunction*>(HEOS.Reducing.get());
1039 return HEOS.rhomolar_reducing() * GERG->dTrdxi__constxj(HEOS.mole_fractions, i, xN_flag)
1040 + HEOS.T_reducing() * GERG->drhormolardxi__constxj(HEOS.mole_fractions, i, xN_flag);
1041}
1042
1043CoolPropDbl MixtureDerivatives::d2_rhorTr_dxidxj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
1044 GERG2008ReducingFunction* GERG = static_cast<GERG2008ReducingFunction*>(HEOS.Reducing.get());
1045 return (HEOS.rhomolar_reducing() * GERG->d2Trdxidxj(HEOS.mole_fractions, i, j, xN_flag)
1046 + GERG->drhormolardxi__constxj(HEOS.mole_fractions, j, xN_flag) * GERG->dTrdxi__constxj(HEOS.mole_fractions, i, xN_flag)
1047 + HEOS.T_reducing() * GERG->d2rhormolardxidxj(HEOS.mole_fractions, i, j, xN_flag)
1048 + GERG->dTrdxi__constxj(HEOS.mole_fractions, j, xN_flag) * GERG->drhormolardxi__constxj(HEOS.mole_fractions, i, xN_flag));
1049}
1050
1051CoolPropDbl MixtureDerivatives::d2psi_dDelta2(HelmholtzEOSMixtureBackend& HEOS, x_N_dependency_flag xN_flag) {
1052 return HEOS.rhomolar_reducing() * HEOS.gas_constant() * HEOS.T() * (HEOS.delta() * d2alpha_dDelta2(HEOS) + 2 * dalpha_dDelta(HEOS));
1053}
1054
1055CoolPropDbl MixtureDerivatives::d2psi_dDelta_dTau(HelmholtzEOSMixtureBackend& HEOS, x_N_dependency_flag xN_flag) {
1056 return HEOS.rhomolar_reducing() * HEOS.gas_constant() * HEOS.T() / HEOS.tau()
1057 * (HEOS.tau() * dalpha_dTau(HEOS) - alpha(HEOS) - HEOS.delta() * dalpha_dDelta(HEOS)
1058 + HEOS.tau() * HEOS.delta() * d2alpha_dDelta_dTau(HEOS));
1059}
1060CoolPropDbl MixtureDerivatives::d2psir_dDelta_dTau(HelmholtzEOSMixtureBackend& HEOS, x_N_dependency_flag xN_flag) {
1061 return HEOS.rhomolar_reducing() * HEOS.gas_constant() * HEOS.T() / HEOS.tau()
1062 * (HEOS.tau() * dalphar_dTau(HEOS) - alphar(HEOS) - HEOS.delta() * dalphar_dDelta(HEOS)
1063 + HEOS.tau() * HEOS.delta() * d2alphar_dDelta_dTau(HEOS));
1064}
1065CoolPropDbl MixtureDerivatives::d2psi_dTau2(HelmholtzEOSMixtureBackend& HEOS, x_N_dependency_flag xN_flag) {
1066 double tau = HEOS.tau();
1067 return HEOS.rhomolar_reducing() * HEOS.delta() * HEOS.gas_constant() * HEOS.T() / POW2(tau)
1068 * (POW2(tau) * d2alpha_dTau2(HEOS) - 2 * tau * dalpha_dTau(HEOS) + 2 * alpha(HEOS));
1069}
1070CoolPropDbl MixtureDerivatives::d2psir_dTau2(HelmholtzEOSMixtureBackend& HEOS, x_N_dependency_flag xN_flag) {
1071 double tau = HEOS.tau();
1072 return HEOS.rhomolar_reducing() * HEOS.delta() * HEOS.gas_constant() * HEOS.T() / POW2(tau)
1073 * (POW2(tau) * d2alphar_dTau2(HEOS) - 2 * tau * dalphar_dTau(HEOS) + 2 * alphar(HEOS));
1074}
1075CoolPropDbl MixtureDerivatives::d2psi_dxi_dDelta(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
1076 return HEOS.gas_constant() / HEOS.tau()
1077 * (d_rhorTr_dxi(HEOS, i, xN_flag) * (HEOS.delta() * dalpha_dDelta(HEOS) + alpha(HEOS))
1078 + HEOS.rhomolar_reducing() * HEOS.T_reducing() * (HEOS.delta() * d2alpha_dxi_dDelta(HEOS, i, xN_flag) + dalpha_dxi(HEOS, i, xN_flag)));
1079}
1080CoolPropDbl MixtureDerivatives::d2psi_dxi_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
1081 return HEOS.delta() * HEOS.gas_constant() / POW2(HEOS.tau())
1082 * (d_rhorTr_dxi(HEOS, i, xN_flag) * (HEOS.tau() * dalpha_dTau(HEOS) - alpha(HEOS))
1083 + HEOS.rhomolar_reducing() * HEOS.T_reducing() * (HEOS.tau() * d2alpha_dxi_dTau(HEOS, i, xN_flag) - dalpha_dxi(HEOS, i, xN_flag)));
1084}
1085CoolPropDbl MixtureDerivatives::d2psir_dxi_dTau(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, x_N_dependency_flag xN_flag) {
1086 return HEOS.delta() * HEOS.gas_constant() / POW2(HEOS.tau())
1087 * (d_rhorTr_dxi(HEOS, i, xN_flag) * (HEOS.tau() * dalphar_dTau(HEOS) - alphar(HEOS))
1088 + HEOS.rhomolar_reducing() * HEOS.T_reducing() * (HEOS.tau() * d2alphar_dxi_dTau(HEOS, i, xN_flag) - dalphar_dxi(HEOS, i, xN_flag)));
1089}
1090CoolPropDbl MixtureDerivatives::d2psi_dxi_dxj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
1091 return HEOS.delta() * HEOS.gas_constant() / HEOS.tau()
1092 * (alpha(HEOS) * d2_rhorTr_dxidxj(HEOS, i, j, xN_flag) + dalpha_dxi(HEOS, i, xN_flag) * d_rhorTr_dxi(HEOS, j, xN_flag)
1093 + dalpha_dxi(HEOS, j, xN_flag) * d_rhorTr_dxi(HEOS, i, xN_flag)
1094 + HEOS.rhomolar_reducing() * HEOS.T_reducing() * d2alphadxidxj(HEOS, i, j, xN_flag));
1095}
1096CoolPropDbl MixtureDerivatives::d2psir_dxi_dxj(HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, x_N_dependency_flag xN_flag) {
1097 return HEOS.delta() * HEOS.gas_constant() / HEOS.tau()
1098 * (alphar(HEOS) * d2_rhorTr_dxidxj(HEOS, i, j, xN_flag) + dalphar_dxi(HEOS, i, xN_flag) * d_rhorTr_dxi(HEOS, j, xN_flag)
1099 + dalphar_dxi(HEOS, j, xN_flag) * d_rhorTr_dxi(HEOS, i, xN_flag)
1100 + HEOS.rhomolar_reducing() * HEOS.T_reducing() * d2alphardxidxj(HEOS, i, j, xN_flag));
1101}
1102
1103} /* namespace CoolProp */
1104
1105#ifdef ENABLE_CATCH
1106# include <catch2/catch_all.hpp>
1107
1108# include "Backends/Cubics/CubicBackend.h"
1109# include "Backends/Cubics/VTPRBackend.h"
1110
1111using namespace CoolProp;
1112
1113double mix_deriv_err_func(double numeric, double analytic) {
1114 if (std::abs(analytic) < DBL_EPSILON) {
1115 return std::abs(numeric - analytic);
1116 } else {
1117 return std::abs(numeric / analytic - 1);
1118 }
1119}
1120
1121typedef CoolProp::MixtureDerivatives MD;
1122
1123enum derivative
1124{
1125 NO_DERIV = 0,
1126 TAU,
1127 DELTA,
1128 XI,
1129 XJ,
1130 XK,
1131 T_CONSTP,
1132 P_CONSTT,
1133 T_CONSTRHO,
1134 RHO_CONSTT,
1135 CONST_TAU_DELTA,
1136 CONST_TRHO
1137};
1138
1139typedef double (*zero_mole_fraction_pointer)(CoolProp::HelmholtzEOSMixtureBackend& HEOS, CoolProp::x_N_dependency_flag xN_flag);
1140typedef double (*one_mole_fraction_pointer)(CoolProp::HelmholtzEOSMixtureBackend& HEOS, std::size_t i, CoolProp::x_N_dependency_flag xN_flag);
1141typedef double (*two_mole_fraction_pointer)(CoolProp::HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j,
1142 CoolProp::x_N_dependency_flag xN_flag);
1143typedef double (*three_mole_fraction_pointer)(CoolProp::HelmholtzEOSMixtureBackend& HEOS, std::size_t i, std::size_t j, std::size_t k,
1144 CoolProp::x_N_dependency_flag xN_flag);
1145
1146template <class backend>
1147class DerivativeFixture
1148{
1149 public:
1150 shared_ptr<CoolProp::HelmholtzEOSMixtureBackend> HEOS, HEOS_plus_tau, HEOS_minus_tau, HEOS_plus_delta, HEOS_minus_delta, HEOS_plus_2tau, HEOS_minus_2tau, HEOS_plus_2delta, HEOS_minus_2delta, HEOS_plus_T__constp,
1151 HEOS_minus_T__constp, HEOS_plus_p__constT, HEOS_minus_p__constT, HEOS_plus_T__constrho, HEOS_minus_T__constrho, HEOS_plus_rho__constT,
1152 HEOS_minus_rho__constT;
1153 std::vector<shared_ptr<CoolProp::HelmholtzEOSMixtureBackend>> HEOS_plus_z, HEOS_minus_z, HEOS_plus_z__constTrho, HEOS_minus_z__constTrho,
1154 HEOS_plus_n, HEOS_minus_n, HEOS_plus_2z, HEOS_minus_2z, HEOS_plus_2z__constTrho, HEOS_minus_2z__constTrho;
1155 CoolProp::x_N_dependency_flag xN;
1156 double dtau, ddelta, dz, dn, tol, dT, drho, dp;
1157 DerivativeFixture() : xN(XN_INDEPENDENT) {
1158 dtau = 1e-6;
1159 ddelta = 1e-6;
1160 dz = 1e-6;
1161 dn = 1e-6;
1162 dT = 1e-3;
1163 drho = 1e-3;
1164 dp = 1;
1165 tol = 5e-6;
1166 std::vector<std::string> names;
1167 names.push_back("n-Pentane");
1168 names.push_back("Ethane");
1169 names.push_back("n-Propane");
1170 names.push_back("n-Butane");
1171 std::vector<CoolPropDbl> mole_fractions;
1172 mole_fractions.push_back(0.1);
1173 mole_fractions.push_back(0.12);
1174 mole_fractions.push_back(0.18);
1175 mole_fractions.push_back(0.6);
1176 HEOS.reset(new backend(names));
1177 HEOS->set_mole_fractions(mole_fractions);
1178 HEOS->specify_phase(CoolProp::iphase_gas);
1179 HEOS->update_DmolarT_direct(300, 300);
1180
1181 init_state(HEOS_plus_tau);
1182 init_state(HEOS_minus_tau);
1183 init_state(HEOS_plus_2tau);
1184 init_state(HEOS_minus_2tau);
1185 init_state(HEOS_plus_delta);
1186 init_state(HEOS_minus_delta);
1187 init_state(HEOS_plus_2delta);
1188 init_state(HEOS_minus_2delta);
1189 init_state(HEOS_plus_T__constp);
1190 init_state(HEOS_minus_T__constp);
1191 init_state(HEOS_plus_p__constT);
1192 init_state(HEOS_minus_p__constT);
1193 init_state(HEOS_plus_T__constrho);
1194 init_state(HEOS_minus_T__constrho);
1195 init_state(HEOS_plus_rho__constT);
1196 init_state(HEOS_minus_rho__constT);
1197 // Constant composition, varying state
1198 HEOS_plus_tau->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS->rhomolar_reducing(), HEOS->T_reducing() / (HEOS->tau() + dtau));
1199 HEOS_minus_tau->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS->rhomolar_reducing(), HEOS->T_reducing() / (HEOS->tau() - dtau));
1200 HEOS_plus_2tau->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS->rhomolar_reducing(), HEOS->T_reducing() / (HEOS->tau() + 2*dtau));
1201 HEOS_minus_2tau->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS->rhomolar_reducing(), HEOS->T_reducing() / (HEOS->tau() - 2*dtau));
1202 HEOS_plus_delta->update(CoolProp::DmolarT_INPUTS, (HEOS->delta() + ddelta) * HEOS->rhomolar_reducing(), HEOS->T_reducing() / HEOS->tau());
1203 HEOS_minus_delta->update(CoolProp::DmolarT_INPUTS, (HEOS->delta() - ddelta) * HEOS->rhomolar_reducing(), HEOS->T_reducing() / HEOS->tau());
1204 HEOS_plus_2delta->update(CoolProp::DmolarT_INPUTS, (HEOS->delta() + 2*ddelta) * HEOS->rhomolar_reducing(), HEOS->T_reducing() / HEOS->tau());
1205 HEOS_minus_2delta->update(CoolProp::DmolarT_INPUTS, (HEOS->delta() - 2*ddelta) * HEOS->rhomolar_reducing(), HEOS->T_reducing() / HEOS->tau());
1206 HEOS_plus_T__constp->update(CoolProp::PT_INPUTS, HEOS->p(), HEOS->T() + dT);
1207 HEOS_minus_T__constp->update(CoolProp::PT_INPUTS, HEOS->p(), HEOS->T() - dT);
1208 HEOS_plus_p__constT->update(CoolProp::PT_INPUTS, HEOS->p() + dp, HEOS->T());
1209 HEOS_minus_p__constT->update(CoolProp::PT_INPUTS, HEOS->p() - dp, HEOS->T());
1210 HEOS_plus_T__constrho->update(CoolProp::DmolarT_INPUTS, HEOS->rhomolar(), HEOS->T() + dT);
1211 HEOS_minus_T__constrho->update(CoolProp::DmolarT_INPUTS, HEOS->rhomolar(), HEOS->T() - dT);
1212 HEOS_plus_rho__constT->update(CoolProp::DmolarT_INPUTS, HEOS->rhomolar() + drho, HEOS->T());
1213 HEOS_minus_rho__constT->update(CoolProp::DmolarT_INPUTS, HEOS->rhomolar() - drho, HEOS->T());
1214
1215 // Varying mole fractions
1216 HEOS_plus_z.resize(4);
1217 HEOS_minus_z.resize(4);
1218 HEOS_plus_z__constTrho.resize(4);
1219 HEOS_minus_z__constTrho.resize(4);
1220 HEOS_plus_2z.resize(4);
1221 HEOS_minus_2z.resize(4);
1222 HEOS_plus_2z__constTrho.resize(4);
1223 HEOS_minus_2z__constTrho.resize(4);
1224 for (int i = 0; i < HEOS_plus_z.size(); ++i) {
1225 init_state(HEOS_plus_z[i]);
1226 init_state(HEOS_plus_2z[i]);
1227 init_state(HEOS_plus_z__constTrho[i]);
1228 std::vector<double> zz = HEOS->get_mole_fractions(), zz2 = HEOS->get_mole_fractions();
1229 zz[i] += dz;
1230 zz2[i] += 2 * dz;
1231 if (xN == CoolProp::XN_DEPENDENT) {
1232 zz[zz.size() - 1] -= dz;
1233 zz2[zz2.size() - 1] -= 2 * dz;
1234 }
1235 HEOS_plus_z[i]->set_mole_fractions(zz);
1236 HEOS_plus_z[i]->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS_plus_z[i]->rhomolar_reducing(),
1237 HEOS_plus_z[i]->T_reducing() / HEOS->tau());
1238 HEOS_plus_2z[i]->set_mole_fractions(zz2);
1239 HEOS_plus_2z[i]->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS_plus_2z[i]->rhomolar_reducing(),
1240 HEOS_plus_2z[i]->T_reducing() / HEOS->tau());
1241 HEOS_plus_z__constTrho[i]->set_mole_fractions(zz);
1242 HEOS_plus_z__constTrho[i]->update(CoolProp::DmolarT_INPUTS, HEOS->rhomolar(), HEOS->T());
1243 }
1244 for (int i = 0; i < HEOS_minus_z.size(); ++i) {
1245 init_state(HEOS_minus_z[i]);
1246 init_state(HEOS_minus_2z[i]);
1247 init_state(HEOS_minus_z__constTrho[i]);
1248 std::vector<double> zz = HEOS->get_mole_fractions(), zz2 = HEOS->get_mole_fractions();
1249 zz[i] -= dz;
1250 zz2[i] -= 2 * dz;
1251 if (xN == CoolProp::XN_DEPENDENT) {
1252 zz[zz.size() - 1] += dz;
1253 zz2[zz2.size() - 1] += 2 * dz;
1254 }
1255 HEOS_minus_z[i]->set_mole_fractions(zz);
1256 HEOS_minus_z[i]->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS_minus_z[i]->rhomolar_reducing(),
1257 HEOS_minus_z[i]->T_reducing() / HEOS->tau());
1258 HEOS_minus_2z[i]->set_mole_fractions(zz2);
1259 HEOS_minus_2z[i]->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS_minus_2z[i]->rhomolar_reducing(),
1260 HEOS_minus_2z[i]->T_reducing() / HEOS->tau());
1261 HEOS_minus_z__constTrho[i]->set_mole_fractions(zz);
1262 HEOS_minus_z__constTrho[i]->update(CoolProp::DmolarT_INPUTS, HEOS->rhomolar(), HEOS->T());
1263 }
1264
1265 // Varying mole numbers
1266 HEOS_plus_n.resize(4);
1267 HEOS_minus_n.resize(4);
1268 for (int i = 0; i < HEOS_plus_n.size(); ++i) {
1269 init_state(HEOS_plus_n[i]);
1270 std::vector<double> zz = HEOS->get_mole_fractions();
1271 zz[i] += dn;
1272 for (int j = 0; j < HEOS_minus_n.size(); ++j) {
1273 zz[i] /= 1 + dn;
1274 }
1275 HEOS_plus_n[i]->set_mole_fractions(zz);
1276 HEOS_plus_n[i]->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS_plus_n[i]->rhomolar_reducing(),
1277 HEOS_plus_n[i]->T_reducing() / HEOS->tau());
1278 }
1279 for (int i = 0; i < HEOS_plus_z.size(); ++i) {
1280 init_state(HEOS_minus_n[i]);
1281 std::vector<double> zz = HEOS->get_mole_fractions();
1282 zz[i] -= dn;
1283 for (int j = 0; j < HEOS_minus_n.size(); ++j) {
1284 zz[i] /= 1 - dn;
1285 }
1286 HEOS_minus_n[i]->set_mole_fractions(zz);
1287 HEOS_minus_n[i]->update(CoolProp::DmolarT_INPUTS, HEOS->delta() * HEOS_minus_n[i]->rhomolar_reducing(),
1288 HEOS_minus_n[i]->T_reducing() / HEOS->tau());
1289 }
1290 };
1291 void init_state(shared_ptr<CoolProp::HelmholtzEOSMixtureBackend>& other) {
1292 other.reset(HEOS->get_copy());
1293 other->set_mole_fractions(HEOS->get_mole_fractions());
1294 other->specify_phase(CoolProp::iphase_gas); // Something homogeneous
1295 }
1296 void zero(const std::string& name, zero_mole_fraction_pointer f, zero_mole_fraction_pointer g = NULL, derivative wrt = NO_DERIV) {
1297 double analytic = f(*HEOS, xN);
1298 double numeric = 0;
1299 if (wrt == TAU) {
1300 numeric = (g(*HEOS_plus_tau, xN) - g(*HEOS_minus_tau, xN)) / (2 * dtau);
1301 } else if (wrt == DELTA) {
1302 numeric = (g(*HEOS_plus_delta, xN) - g(*HEOS_minus_delta, xN)) / (2 * ddelta);
1303 }
1304 CAPTURE(name);
1305 CAPTURE(analytic);
1306 CAPTURE(numeric);
1307 CAPTURE(xN);
1308 double error = mix_deriv_err_func(numeric, analytic);
1309 CAPTURE(error);
1310 CHECK(error < tol);
1311 }
1312 void one(const std::string& name, one_mole_fraction_pointer f, one_mole_fraction_pointer g = NULL, derivative wrt = NO_DERIV) {
1313 for (int i = 0; i < 4; ++i) {
1314 double analytic = f(*HEOS, i, xN);
1315 double numeric = 0;
1316 if (wrt == TAU) {
1317 numeric = (g(*HEOS_plus_tau, i, xN) - g(*HEOS_minus_tau, i, xN)) / (2 * dtau);
1318 } else if (wrt == DELTA) {
1319 numeric = (g(*HEOS_plus_delta, i, xN) - g(*HEOS_minus_delta, i, xN)) / (2 * ddelta);
1320 } else if (wrt == T_CONSTP) {
1321 numeric = (g(*HEOS_plus_T__constp, i, xN) - g(*HEOS_minus_T__constp, i, xN)) / (2 * dT);
1322 } else if (wrt == P_CONSTT) {
1323 numeric = (g(*HEOS_plus_p__constT, i, xN) - g(*HEOS_minus_p__constT, i, xN)) / (2 * dp);
1324 } else if (wrt == T_CONSTRHO) {
1325 double g1 = g(*HEOS_plus_T__constrho, i, xN), g2 = g(*HEOS_minus_T__constrho, i, xN);
1326 numeric = (g1 - g2) / (2 * dT);
1327 } else if (wrt == RHO_CONSTT) {
1328 numeric = (g(*HEOS_plus_rho__constT, i, xN) - g(*HEOS_minus_rho__constT, i, xN)) / (2 * drho);
1329 }
1330 CAPTURE(name);
1331 CAPTURE(i);
1332 CAPTURE(analytic);
1333 CAPTURE(numeric);
1334 CAPTURE(xN);
1335 double error = mix_deriv_err_func(numeric, analytic);
1336 CAPTURE(error);
1337 CHECK(error < tol);
1338 }
1339 };
1340 void one_comp(const std::string& name, one_mole_fraction_pointer f, zero_mole_fraction_pointer g, derivative wrt = CONST_TAU_DELTA) {
1341 for (int i = 0; i < 4; ++i) {
1342 double analytic;
1343 CHECK_NOTHROW(analytic = f(*HEOS, i, xN));
1344 double numeric = -10000;
1345 if (wrt == CONST_TAU_DELTA) {
1346 if (HEOS->get_mole_fractions()[i] > dz) {
1347 CHECK_NOTHROW(numeric = (g(*HEOS_plus_z[i], xN) - g(*HEOS_minus_z[i], xN)) / (2 * dz));
1348 } else {
1349 CHECK_NOTHROW(numeric = (-3 * g(*HEOS, xN) + 4 * g(*HEOS_plus_z[i], xN) - g(*HEOS_plus_2z[i], xN)) / (2 * dz));
1350 }
1351 } else if (wrt == CONST_TRHO) {
1352 CHECK_NOTHROW(numeric = (g(*HEOS_plus_z__constTrho[i], xN) - g(*HEOS_minus_z__constTrho[i], xN)) / (2 * dz));
1353 }
1354
1355 CAPTURE(name);
1356 CAPTURE(i);
1357 CAPTURE(analytic);
1358 CAPTURE(numeric);
1359 CAPTURE(xN);
1360 double error = mix_deriv_err_func(numeric, analytic);
1361 CAPTURE(error);
1362 CHECK(error < tol);
1363 }
1364 }
1365 void two(const std::string& name, two_mole_fraction_pointer f, two_mole_fraction_pointer g = NULL, derivative wrt = NO_DERIV) {
1366 for (int i = 0; i < 4; ++i) {
1367 for (int j = 0; j < 4; ++j) {
1368 double analytic = f(*HEOS, i, j, xN);
1369 bool is_other = false;
1370 double numeric = 0;
1371 if (wrt == TAU) {
1372 numeric = (g(*HEOS_plus_tau, i, j, xN) - g(*HEOS_minus_tau, i, j, xN)) / (2 * dtau);
1373 is_other = true;
1374 } else if (wrt == DELTA) {
1375 numeric = (g(*HEOS_plus_delta, i, j, xN) - g(*HEOS_minus_delta, i, j, xN)) / (2 * ddelta);
1376 is_other = true;
1377 }
1378 CAPTURE(name);
1379 CAPTURE(i);
1380 CAPTURE(j);
1381 CAPTURE(analytic);
1382 CAPTURE(numeric);
1383 CAPTURE(xN);
1384 double error = mix_deriv_err_func(numeric, analytic);
1385 CAPTURE(error);
1386 CHECK(error < tol);
1387 }
1388 }
1389 }
1390 void two_comp(const std::string& name, two_mole_fraction_pointer f, one_mole_fraction_pointer g, derivative wrt = NO_DERIV) {
1391 for (int i = 0; i < 4; ++i) {
1392 for (int j = 0; j < 4; ++j) {
1393 double analytic = f(*HEOS, i, j, xN);
1394 double numeric = 500;
1395 if (HEOS->get_mole_fractions()[j] > 2 * dz) {
1396 // Second order centered difference in composition
1397 CHECK_NOTHROW(numeric = (g(*HEOS_plus_z[j], i, xN) - g(*HEOS_minus_z[j], i, xN)) / (2 * dz));
1398 } else {
1399 // Forward difference in composition
1400 CHECK_NOTHROW(numeric = (-3 * g(*HEOS, i, xN) + 4 * g(*HEOS_plus_z[j], i, xN) - g(*HEOS_plus_2z[j], i, xN)) / (2 * dz));
1401 }
1402 CAPTURE(name);
1403 CAPTURE(i);
1404 CAPTURE(j);
1405 CAPTURE(analytic);
1406 CAPTURE(numeric);
1407 CAPTURE(xN);
1408 double error = mix_deriv_err_func(numeric, analytic);
1409 CAPTURE(error);
1410 CHECK(error < tol);
1411 }
1412 }
1413 }
1414 void three(const std::string& name, three_mole_fraction_pointer f, three_mole_fraction_pointer g = NULL, derivative wrt = NO_DERIV) {
1415 for (int i = 0; i < 4; ++i) {
1416 for (int j = 0; j < 4; ++j) {
1417 for (int k = 0; k < 4; ++k) {
1418 double analytic = f(*HEOS, i, j, k, xN);
1419 double numeric = 0;
1420 if (wrt == TAU) {
1421 numeric = (-1.0/12.0*g(*HEOS_plus_2tau, i, j, k, xN) + 2.0/3.0*g(*HEOS_plus_tau, i, j, k, xN) -2.0/3.0*g(*HEOS_minus_tau, i, j, k, xN) +1.0/12.0*g(*HEOS_minus_2tau, i, j, k, xN)) / dtau;
1422 } else if (wrt == DELTA) {
1423 numeric = (-1.0/12.0*g(*HEOS_plus_2delta, i, j, k, xN) + 2.0/3.0*g(*HEOS_plus_delta, i, j, k, xN) -2.0/3.0*g(*HEOS_minus_delta, i, j, k, xN) +1.0/12.0*g(*HEOS_minus_2delta, i, j, k, xN)) / ddelta;
1424 }
1425 CAPTURE(name);
1426 CAPTURE(i);
1427 CAPTURE(j);
1428 CAPTURE(analytic);
1429 CAPTURE(numeric);
1430 CAPTURE(xN);
1431 double error = mix_deriv_err_func(numeric, analytic);
1432 CAPTURE(error);
1433 CHECK(error < tol);
1434 }
1435 }
1436 }
1437 }
1438 void three_comp(const std::string& name, three_mole_fraction_pointer f, two_mole_fraction_pointer g, derivative wrt = NO_DERIV) {
1439 for (int i = 0; i < 4; ++i) {
1440 for (int j = 0; j < 4; ++j) {
1441 for (int k = 0; k < 4; ++k) {
1442 double analytic = f(*HEOS, i, j, k, xN);
1443 double numeric;
1444 if (HEOS->get_mole_fractions()[i] > 2 * dz) {
1445 CHECK_NOTHROW(numeric = (g(*HEOS_plus_z[k], i, j, xN) - g(*HEOS_minus_z[k], i, j, xN)) / (2 * dz));
1446 } else {
1447 CHECK_NOTHROW(numeric =
1448 (-3 * g(*HEOS, i, j, xN) + 4 * g(*HEOS_plus_z[k], i, j, xN) - g(*HEOS_plus_2z[k], i, j, xN)) / (2 * dz));
1449 }
1450 CAPTURE(name);
1451 CAPTURE(i);
1452 CAPTURE(j);
1453 CAPTURE(k);
1454 CAPTURE(analytic);
1455 CAPTURE(numeric);
1456 CAPTURE(xN);
1457 double error = mix_deriv_err_func(numeric, analytic);
1458 CAPTURE(error);
1459 CHECK(error < tol);
1460 }
1461 }
1462 }
1463 }
1464 Eigen::MatrixXd get_matrix(CoolProp::HelmholtzEOSMixtureBackend& HEOS, const std::string name) {
1465 Eigen::MatrixXd Lstar = MixtureDerivatives::Lstar(HEOS, xN);
1466 if (name == "Lstar") {
1467 return Lstar;
1468 } else if (name == "Mstar") {
1469 return MixtureDerivatives::Mstar(HEOS, xN, Lstar);
1470 } else {
1471 throw ValueError("Must be one of Lstar or Mstar");
1472 }
1473 }
1474 void matrix_derivative(const std::string name, const std::string& wrt) {
1475 CAPTURE(name);
1476 CAPTURE(wrt);
1477 double err = 10000;
1478 Eigen::MatrixXd analytic, numeric, Lstar, dLstar_dTau, dLstar_dDelta;
1479 if (name == "Mstar") {
1480 Lstar = MixtureDerivatives::Lstar(*HEOS, xN);
1481 dLstar_dDelta = MixtureDerivatives::dLstar_dX(*HEOS, xN, CoolProp::iDelta);
1482 dLstar_dTau = MixtureDerivatives::dLstar_dX(*HEOS, xN, CoolProp::iTau);
1483 }
1484 if (wrt == "tau") {
1485 if (name == "Lstar") {
1486 analytic = MixtureDerivatives::dLstar_dX(*HEOS, xN, CoolProp::iTau);
1487 } else if (name == "Mstar") {
1488 analytic = MixtureDerivatives::dMstar_dX(*HEOS, xN, CoolProp::iTau, Lstar, dLstar_dTau);
1489 } else {
1490 throw ValueError("argument not understood");
1491 }
1492 numeric = (get_matrix(*HEOS_plus_tau, name) - get_matrix(*HEOS_minus_tau, name)) / (2 * dtau);
1493 } else if (wrt == "delta") {
1494 if (name == "Lstar") {
1495 analytic = MixtureDerivatives::dLstar_dX(*HEOS, xN, CoolProp::iDelta);
1496 } else if (name == "Mstar") {
1497 analytic = MixtureDerivatives::dMstar_dX(*HEOS, xN, CoolProp::iDelta, Lstar, dLstar_dDelta);
1498 } else {
1499 throw ValueError("argument not understood");
1500 }
1501 numeric = (get_matrix(*HEOS_plus_delta, name) - get_matrix(*HEOS_minus_delta, name)) / (2 * ddelta);
1502 } else {
1503 throw ValueError("argument not understood");
1504 }
1505 CAPTURE(analytic);
1506 CAPTURE(numeric);
1507 Eigen::MatrixXd rel_error = ((analytic - numeric).cwiseQuotient(analytic));
1508 CAPTURE(rel_error);
1509 err = rel_error.squaredNorm();
1510 CAPTURE(err);
1511 CHECK(err < 1e-8);
1512 }
1513
1514 void run_checks() {
1515
1516 two_comp("d_PSI_rho_dxj", MD::d_PSI_rho_dxj, MD::PSI_rho);
1517 two_comp("d_PSI_T_dxj", MD::d_PSI_T_dxj, MD::PSI_T);
1518
1519 one("d_ndalphardni_dDelta", MD::d_ndalphardni_dDelta, MD::ndalphar_dni__constT_V_nj, DELTA);
1520 one("d2_ndalphardni_dDelta2", MD::d2_ndalphardni_dDelta2, MD::d_ndalphardni_dDelta, DELTA);
1521 one("d3_ndalphardni_dDelta3", MD::d3_ndalphardni_dDelta3, MD::d2_ndalphardni_dDelta2, DELTA);
1522 one("d_ndalphardni_dTau", MD::d_ndalphardni_dTau, MD::ndalphar_dni__constT_V_nj, TAU);
1523 one("d2_ndalphardni_dTau2", MD::d2_ndalphardni_dTau2, MD::d_ndalphardni_dTau, TAU);
1524 one("d3_ndalphardni_dTau3", MD::d3_ndalphardni_dTau3, MD::d2_ndalphardni_dTau2, TAU);
1525 one("d2_ndalphardni_dDelta_dTau", MD::d2_ndalphardni_dDelta_dTau, MD::d_ndalphardni_dDelta, TAU);
1526 one("d3_ndalphardni_dDelta2_dTau", MD::d3_ndalphardni_dDelta2_dTau, MD::d2_ndalphardni_dDelta2, TAU);
1527 one("d3_ndalphardni_dDelta_dTau2", MD::d3_ndalphardni_dDelta_dTau2, MD::d2_ndalphardni_dDelta_dTau, TAU);
1528
1529 zero("dalphar_dDelta", MD::dalphar_dDelta, MD::alphar, DELTA);
1530 zero("d2alphar_dDelta2", MD::d2alphar_dDelta2, MD::dalphar_dDelta, DELTA);
1531 zero("dalphar_dTau", MD::dalphar_dTau, MD::alphar, TAU);
1532 zero("d2alphar_dTau2", MD::d2alphar_dTau2, MD::dalphar_dTau, TAU);
1533
1534 zero("dalpha0_dDelta", MD::dalpha0_dDelta, MD::alpha0, DELTA);
1535 zero("d2alpha0_dDelta2", MD::d2alpha0_dDelta2, MD::dalpha0_dDelta, DELTA);
1536 zero("dalpha0_dTau", MD::dalpha0_dTau, MD::alpha0, TAU);
1537 zero("d2alpha0_dTau2", MD::d2alpha0_dTau2, MD::dalpha0_dTau, TAU);
1538
1539 one_comp("dalpha0_dxi", MD::dalpha0_dxi, MD::alpha0);
1540 one("d2alpha0_dxi_dDelta", MD::d2alpha0_dxi_dDelta, MD::dalpha0_dxi, DELTA);
1541 one("d2alpha0_dxi_dTau", MD::d2alpha0_dxi_dTau, MD::dalpha0_dxi, TAU);
1542 two_comp("d2alpha0dxidxj", MD::d2alpha0dxidxj, MD::dalpha0_dxi);
1543
1544 one_comp("dalpha_dxi", MD::dalpha_dxi, MD::alpha);
1545 one("d2alpha_dxi_dDelta", MD::d2alpha_dxi_dDelta, MD::dalpha_dxi, DELTA);
1546 one("d2alpha_dxi_dTau", MD::d2alpha_dxi_dTau, MD::dalpha_dxi, TAU);
1547 two_comp("d2alphadxidxj", MD::d2alphadxidxj, MD::dalpha_dxi);
1548
1549 zero("dpsi_dDelta", MD::dpsi_dDelta, MD::psi, DELTA);
1550 zero("dpsi_dTau", MD::dpsi_dTau, MD::psi, TAU);
1551 zero("d2psi_dDelta2", MD::d2psi_dDelta2, MD::dpsi_dDelta, DELTA);
1552 zero("d2psi_dDelta_dTau", MD::d2psi_dDelta_dTau, MD::dpsi_dDelta, TAU);
1553 zero("d2psi_dTau2", MD::d2psi_dTau2, MD::dpsi_dTau, TAU);
1554 one_comp("dpsi_dxi", MD::dpsi_dxi, MD::psi);
1555 one_comp("d2psi_dxi_dDelta", MD::d2psi_dxi_dDelta, MD::dpsi_dDelta);
1556 one_comp("d2psi_dxi_dTau", MD::d2psi_dxi_dTau, MD::dpsi_dTau);
1557 two_comp("d2psi_dxi_dxj", MD::d2psi_dxi_dxj, MD::dpsi_dxi);
1558
1559 //two_comp("d_ndalphardni_dxj__constT_V_xi", MD::d_ndalphardni_dxj__constT_V_xi, MD::ndalphar_dni__constT_V_nj);
1560
1561 one_comp("dalphar_dxi", MD::dalphar_dxi, MD::alphar);
1562 two_comp("d2alphardxidxj", MD::d2alphardxidxj, MD::dalphar_dxi);
1563 three_comp("d3alphardxidxjdxk", MD::d3alphardxidxjdxk, MD::d2alphardxidxj);
1564 one("d2alphar_dxi_dTau", MD::d2alphar_dxi_dTau, MD::dalphar_dxi, TAU);
1565 one("d2alphar_dxi_dDelta", MD::d2alphar_dxi_dDelta, MD::dalphar_dxi, DELTA);
1566 one("d3alphar_dxi_dDelta2", MD::d3alphar_dxi_dDelta2, MD::d2alphar_dxi_dDelta, DELTA);
1567 one("d3alphar_dxi_dTau2", MD::d3alphar_dxi_dTau2, MD::d2alphar_dxi_dTau, TAU);
1568 one("d4alphar_dxi_dTau3", MD::d4alphar_dxi_dTau3, MD::d3alphar_dxi_dTau2, TAU);
1569 one("d3alphar_dxi_dDelta_dTau", MD::d3alphar_dxi_dDelta_dTau, MD::d2alphar_dxi_dDelta, TAU);
1570 one("d4alphar_dxi_dDelta_dTau2", MD::d4alphar_dxi_dDelta_dTau2, MD::d3alphar_dxi_dDelta_dTau, TAU);
1571 one("d4alphar_dxi_dDelta2_dTau", MD::d4alphar_dxi_dDelta2_dTau, MD::d3alphar_dxi_dDelta2, TAU);
1572 two("d3alphar_dxi_dxj_dDelta", MD::d3alphar_dxi_dxj_dDelta, MD::d2alphardxidxj, DELTA);
1573 two("d4alphar_dxi_dxj_dDelta2", MD::d4alphar_dxi_dxj_dDelta2, MD::d3alphar_dxi_dxj_dDelta, DELTA);
1574 two("d4alphar_dxi_dxj_dDelta_dTau", MD::d4alphar_dxi_dxj_dDelta_dTau, MD::d3alphar_dxi_dxj_dDelta, TAU);
1575 two("d3alphar_dxi_dxj_dTau", MD::d3alphar_dxi_dxj_dTau, MD::d2alphardxidxj, TAU);
1576 two("d4alphar_dxi_dxj_dTau2", MD::d4alphar_dxi_dxj_dTau2, MD::d3alphar_dxi_dxj_dTau, TAU);
1577 one_comp("d_dalpharddelta_dxj__constT_V_xi", MD::d_dalpharddelta_dxj__constT_V_xi, MD::dalphar_dDelta, CONST_TRHO);
1578
1579 two_comp("d_ndalphardni_dxj__constdelta_tau_xi", MD::d_ndalphardni_dxj__constdelta_tau_xi, MD::ndalphar_dni__constT_V_nj);
1580 two("d2_ndalphardni_dxj_dDelta__consttau_xi", MD::d2_ndalphardni_dxj_dDelta__consttau_xi, MD::d_ndalphardni_dxj__constdelta_tau_xi, DELTA);
1581 two("d3_ndalphardni_dxj_dDelta2__consttau_xi", MD::d3_ndalphardni_dxj_dDelta2__consttau_xi, MD::d2_ndalphardni_dxj_dDelta__consttau_xi,
1582 DELTA);
1583 two("d2_ndalphardni_dxj_dTau__constdelta_xi", MD::d2_ndalphardni_dxj_dTau__constdelta_xi, MD::d_ndalphardni_dxj__constdelta_tau_xi, TAU);
1584 two("d3_ndalphardni_dxj_dTau2__constdelta_xi", MD::d3_ndalphardni_dxj_dTau2__constdelta_xi, MD::d2_ndalphardni_dxj_dTau__constdelta_xi, TAU);
1585 two("d3_ndalphardni_dxj_dDelta_dTau__constxi", MD::d3_ndalphardni_dxj_dDelta_dTau__constxi, MD::d2_ndalphardni_dxj_dDelta__consttau_xi, TAU);
1586 three_comp("d2_ndalphardni_dxj_dxk__constdelta_tau_xi", MD::d2_ndalphardni_dxj_dxk__constdelta_tau_xi,
1587 MD::d_ndalphardni_dxj__constdelta_tau_xi);
1588 three("d3_ndalphardni_dxj_dxk_dTau__constdelta_xi", MD::d3_ndalphardni_dxj_dxk_dTau__constdelta_xi,
1589 MD::d2_ndalphardni_dxj_dxk__constdelta_tau_xi, TAU);
1590 three("d3_ndalphardni_dxj_dxk_dDelta__consttau_xi", MD::d3_ndalphardni_dxj_dxk_dDelta__consttau_xi,
1591 MD::d2_ndalphardni_dxj_dxk__constdelta_tau_xi, DELTA);
1592
1593 // two("nd_ndalphardni_dnj__constT_V", MD::nd_ndalphardni_dnj__constT_V);
1594 two("d_nd_ndalphardni_dnj_dTau__constdelta_x", MD::d_nd_ndalphardni_dnj_dTau__constdelta_x, MD::nd_ndalphardni_dnj__constT_V, TAU);
1595 two("d2_nd_ndalphardni_dnj_dTau2__constdelta_x", MD::d2_nd_ndalphardni_dnj_dTau2__constdelta_x, MD::d_nd_ndalphardni_dnj_dTau__constdelta_x,
1596 TAU);
1597 two("d_nd_ndalphardni_dnj_dDelta__consttau_x", MD::d_nd_ndalphardni_dnj_dDelta__consttau_x, MD::nd_ndalphardni_dnj__constT_V, DELTA);
1598 two("d2_nd_ndalphardni_dnj_dDelta_dTau__constx", MD::d2_nd_ndalphardni_dnj_dDelta_dTau__constx, MD::d_nd_ndalphardni_dnj_dDelta__consttau_x,
1599 TAU);
1600 two("d2_nd_ndalphardni_dnj_dDelta2__consttau_x", MD::d2_nd_ndalphardni_dnj_dDelta2__consttau_x, MD::d_nd_ndalphardni_dnj_dDelta__consttau_x,
1601 DELTA);
1602 three_comp("d_nd_ndalphardni_dnj_dxk__consttau_delta", MD::d_nd_ndalphardni_dnj_dxk__consttau_delta, MD::nd_ndalphardni_dnj__constT_V);
1603 three("d2_nd_ndalphardni_dnj_dxk_dTau__constdelta", MD::d2_nd_ndalphardni_dnj_dxk_dTau__constdelta,
1604 MD::d_nd_ndalphardni_dnj_dxk__consttau_delta, TAU);
1605 three("d2_nd_ndalphardni_dnj_dxk_dDelta__consttau", MD::d2_nd_ndalphardni_dnj_dxk_dDelta__consttau,
1606 MD::d_nd_ndalphardni_dnj_dxk__consttau_delta, DELTA);
1607
1608 // Xn-dep only two_comp("dln_fugacity_dxj__constT_rho_xi", MD::dln_fugacity_dxj__constT_rho_xi, MD::ln_fugacity);
1609 three("d2_ndln_fugacity_i_dnj_dxk_dDelta__consttau", MD::d2_ndln_fugacity_i_dnj_dxk_dDelta__consttau,
1610 MD::d_ndln_fugacity_i_dnj_ddxk__consttau_delta, DELTA);
1611 three("d2_ndln_fugacity_i_dnj_dxk_dTau__constdelta", MD::d2_ndln_fugacity_i_dnj_dxk_dTau__constdelta,
1612 MD::d_ndln_fugacity_i_dnj_ddxk__consttau_delta, TAU);
1613 two("d2_ndln_fugacity_i_dnj_ddelta_dtau__constx", MD::d2_ndln_fugacity_i_dnj_ddelta_dtau__constx,
1614 MD::d_ndln_fugacity_i_dnj_ddelta__consttau_x, TAU);
1615 two("d_ndln_fugacity_i_dnj_ddelta__consttau_x", MD::d_ndln_fugacity_i_dnj_ddelta__consttau_x, MD::ndln_fugacity_i_dnj__constT_V_xi, DELTA);
1616 two("d_ndln_fugacity_i_dnj_dtau__constdelta_x", MD::d_ndln_fugacity_i_dnj_dtau__constdelta_x, MD::ndln_fugacity_i_dnj__constT_V_xi, TAU);
1617 three_comp("d_ndln_fugacity_i_dnj_ddxk__consttau_delta", MD::d_ndln_fugacity_i_dnj_ddxk__consttau_delta, MD::ndln_fugacity_i_dnj__constT_V_xi,
1618 TAU);
1619 one("dln_fugacity_i_dT__constrho_n", MD::dln_fugacity_i_dT__constrho_n, MD::ln_fugacity, T_CONSTRHO);
1620 one("dln_fugacity_i_drho__constT_n", MD::dln_fugacity_i_drho__constT_n, MD::ln_fugacity, RHO_CONSTT);
1621 one("dln_fugacity_i_dT__constp_n", MD::dln_fugacity_i_dT__constp_n, MD::ln_fugacity, T_CONSTP);
1622 one("dln_fugacity_i_dp__constT_n", MD::dln_fugacity_i_dp__constT_n, MD::ln_fugacity, P_CONSTT);
1623 one("dln_fugacity_coefficient_dT__constp_n", MD::dln_fugacity_coefficient_dT__constp_n, MD::ln_fugacity_coefficient, T_CONSTP);
1624 one("dln_fugacity_coefficient_dp__constT_n", MD::dln_fugacity_coefficient_dp__constT_n, MD::ln_fugacity_coefficient, P_CONSTT);
1625
1626 three_comp("d2_PSI_T_dxj_dxk", MD::d2_PSI_T_dxj_dxk, MD::d_PSI_T_dxj);
1627 three_comp("d2_PSI_rho_dxj_dxk", MD::d2_PSI_rho_dxj_dxk, MD::d_PSI_rho_dxj);
1628
1629 three("d_n2Aijk_dTau", MD::d_n2Aijk_dTau, MD::n2Aijk, TAU);
1630 three("d_n2Aijk_dDelta", MD::d_n2Aijk_dDelta, MD::n2Aijk, DELTA);
1631 two("d_nAij_dTau", MD::d_nAij_dTau, MD::nAij, TAU);
1632 two("d_nAij_dDelta", MD::d_nAij_dDelta, MD::nAij, DELTA);
1633
1634 two_comp("d_nddeltadni_dxj__constdelta_tau", MD::d_nddeltadni_dxj__constdelta_tau, MD::nddeltadni__constT_V_nj);
1635 two_comp("d_ndtaudni_dxj__constdelta_tau", MD::d_ndtaudni_dxj__constdelta_tau, MD::ndtaudni__constT_V_nj);
1636 two("d2_ndtaudni_dxj_dTau__constdelta", MD::d2_ndtaudni_dxj_dTau__constdelta, MD::d_ndtaudni_dxj__constdelta_tau, TAU);
1637
1638 one_comp("dTrdxi__constxj", MD::dTrdxi__constxj, MD::Tr);
1639 one_comp("d2Tr_dxidbetaT", MD::d2Tr_dxidbetaT, MD::dTr_dbetaT);
1640 one_comp("d2Tr_dxidgammaT", MD::d2Tr_dxidgammaT, MD::dTr_dgammaT);
1641 // (??)two_comp("d2Trdxi2__constxj", MD::d2Trdxi2__constxj, MD::dTrdxi__constxj);
1642 two_comp("d2Trdxidxj", MD::d2Trdxidxj, MD::dTrdxi__constxj);
1643 three_comp("d3Trdxidxjdxk", MD::d3Trdxidxjdxk, MD::d2Trdxidxj);
1644 one_comp("dtaudxj__constT_V_xi", MD::dtau_dxj__constT_V_xi, MD::tau, CONST_TRHO);
1645 two_comp("d_ndTrdni_dxj__constxi", MD::d_ndTrdni_dxj__constxi, MD::ndTrdni__constnj);
1646 three_comp("d2_ndTrdni_dxj_dxk__constxi", MD::d2_ndTrdni_dxj_dxk__constxi, MD::d_ndTrdni_dxj__constxi);
1647
1648 one_comp("drhormolardxi__constxj", MD::drhormolardxi__constxj, MD::rhormolar);
1649 one_comp("d2rhormolar_dxidbetaV", MD::d2rhormolar_dxidbetaV, MD::drhormolar_dbetaV);
1650 one_comp("d2rhormolar_dxidgammaV", MD::d2rhormolar_dxidgammaV, MD::drhormolar_dgammaV);
1651 // (??) two_comp("d2rhormolardxi2__constxj", MD::d2rhormolardxi2__constxj, MD::drhormolardxi__constxj);
1652 two_comp("d2rhormolardxidxj", MD::d2rhormolardxidxj, MD::drhormolardxi__constxj);
1653 three_comp("d3rhormolardxidxjdxk", MD::d3rhormolardxidxjdxk, MD::d2rhormolardxidxj);
1654 one_comp("ddelta_dxj__constT_V_xi", MD::ddelta_dxj__constT_V_xi, MD::delta, CONST_TRHO);
1655 two_comp("d_ndrhorbardni_dxj__constxi", MD::d_ndrhorbardni_dxj__constxi, MD::ndrhorbardni__constnj);
1656 three_comp("d2_ndrhorbardni_dxj_dxk__constxi", MD::d2_ndrhorbardni_dxj_dxk__constxi, MD::d_ndrhorbardni_dxj__constxi);
1657
1658 one_comp("dpdxj__constT_V_xi", MD::dpdxj__constT_V_xi, MD::p, CONST_TRHO);
1659
1660 matrix_derivative("Lstar", "tau");
1661 matrix_derivative("Lstar", "delta");
1662 matrix_derivative("Mstar", "tau");
1663 matrix_derivative("Mstar", "delta");
1664 }
1665};
1666
1667TEST_CASE_METHOD(DerivativeFixture<HelmholtzEOSMixtureBackend>, "Check derivatives for HEOS", "[mixture_derivs2]") {
1668 tol = 1e-4; // Relax the tolerance a bit
1669 run_checks();
1670};
1671TEST_CASE_METHOD(DerivativeFixture<PengRobinsonBackend>, "Check derivatives for Peng-Robinson", "[mixture_derivs2]") {
1672 tol = 1e-3; // Relax the tolerance a bit
1673 run_checks();
1674};
1675TEST_CASE_METHOD(DerivativeFixture<SRKBackend>, "Check derivatives for SRK", "[mixture_derivs2]") {
1676 tol = 1e-3; // Relax the tolerance a bit
1677 run_checks();
1678};
1679 // Make sure you set the VTPR UNIFAC path with something like set_config_string(VTPR_UNIFAC_PATH, "/Users/ian/Code/CUBAC/dev/unifaq/");
1680 //TEST_CASE_METHOD(DerivativeFixture<VTPRBackend>, "Check derivatives for VTPR", "[mixture_derivs2]")
1681 //{
1682 // tol = 1e-4; // Relax the tolerance a bit
1683 // run_checks();
1684 //};
1685
1686#endif
Generated on Sat Aug 9 2025 18:35:49 for CoolProp by doxygen 1.9.4