9 for (std::size_t k = 0; k <
N; k++) {
21 for (std::size_t k = 0; k <
N - 1; k++) {
34 for (std::size_t m = 0; m <
N; m++) {
46 for (std::size_t m = 0; m <
N - 1; m++) {
58 for (std::size_t k = 0; k <
N; k++) {
66 for (std::size_t m = 0; m <
N; m++) {
74 for (std::size_t j = 0; j <
N; j++) {
83 for (std::size_t k = 0; k <
N - 1; ++k) {
95 for (std::size_t j = 0; j <
N; j++) {
104 for (std::size_t k = 0; k <
N - 1; ++k) {
125 * (1 -
PSI_rho(x, i, xN_flag));
126 return -1 /
rhormolar(x) * (line1 + line2 + line3 + line4);
141 return 1 /
Tr(x) * (line1 + line2 + line3 + line4);
229 return 2 * pow(rhor, 3) * dvrbardxi * dvrbardxj - pow(rhor, 2) * this->
d2vrmolardxidxj(x, i, j, xN_flag);
235 double line2 = 2 * pow(rhor, 3)
241 return line1 + line2 + line3;
245 const std::vector<CoolPropDbl>& Yc)
const {
247 for (std::size_t i = 0; i <
N; i++) {
249 Yr += xi * xi * Yc[i];
256 for (std::size_t j = i + 1; j <
N; j++) {
264 const STLMatrix& Y_c_ij,
const std::vector<CoolPropDbl>& Yc)
const {
266 for (std::size_t i = 0; i <
N; i++) {
271 for (std::size_t j = i + 1; j <
N; j++) {
278 const STLMatrix& Y_c_ij,
const std::vector<CoolPropDbl>& Yc)
const {
280 for (std::size_t i = 0; i <
N; i++) {
286 for (std::size_t j = i + 1; j <
N; j++) {
287 double xj = x[j], xi = x[i], beta_Y = beta[i][j], beta_Y_squared = beta_Y * beta_Y;
291 double dfYij_dbeta = xi * xj * (-(xi + xj) * (2 * beta_Y * xi)) / pow(beta_Y_squared * xi + xj, 2);
292 dYr_dbeta +=
c_Y_ij(i, j, beta, gamma, Y_c_ij) * dfYij_dbeta +
f_Y_ij(x, i, j, beta) * 2 * gamma[i][j] * Y_c_ij[i][j];
299 const STLMatrix& Y_c_ij,
const std::vector<CoolPropDbl>& Yc,
305 for (std::size_t k = 0; k < i; k++) {
308 for (std::size_t k = i + 1; k <
N; k++) {
317 CoolPropDbl dYr_dxi = 2 * x[i] * Yc[i] - 2 * x[
N - 1] * Yc[
N - 1];
318 for (std::size_t k = 0; k < i; k++) {
321 for (std::size_t k = i + 1; k <
N - 1; k++) {
324 double beta_Y_iN = beta[i][
N - 1], xN = x[
N - 1];
325 dYr_dxi +=
c_Y_ij(i,
N - 1, beta, gamma, Y_c_ij)
326 * (xN * (x[i] + xN) / (pow(beta_Y_iN, 2) * x[i] + xN)
327 + (1 - beta_Y_iN * beta_Y_iN) * x[i] * xN * xN /
POW2(beta_Y_iN * beta_Y_iN * x[i] + xN));
328 for (std::size_t k = 0; k <
N - 1; ++k) {
329 double beta_Y_kN = beta[k][
N - 1], xk = x[k], beta_Y_kN_squared = beta_Y_kN * beta_Y_kN;
331 c_Y_ij(k,
N - 1, beta, gamma, Y_c_ij)
332 * (-xk * (xk + xN) / (beta_Y_kN_squared * xk + xN) + (1 - beta_Y_kN_squared) * xN * xk * xk /
POW2(beta_Y_kN_squared * xk + xN));
345 for (std::size_t k = 0; k < i; k++) {
352 double xk = x[k], beta_Y = beta[k][i], beta_Y_squared = beta_Y * beta_Y;
353 double d2fYkidxidbeta = 2 * beta_Y * pow(xk, 2) * (xi * (xi + xk * (-beta_Y_squared + 1) + xk) - (xi + xk) * (beta_Y_squared * xk + xi))
354 / pow(beta_Y_squared * xk + xi, 3);
355 deriv +=
c_Y_ij(k, i, beta, gamma, Y_c_ij) * d2fYkidxidbeta +
dfYkidxi__constxk(x, k, i, beta) * 2 * gamma[k][i] * Y_c_ij[k][i];
357 for (std::size_t k = i + 1; k <
N; k++) {
363 double xk = x[k], beta_Y = beta[i][k], beta_Y_squared = beta_Y * beta_Y;
364 double d2fYikdxidbeta =
366 * (xi * (-beta_Y_squared * xi + beta_Y_squared * (xi + xk) - xk) - xk * (xi + xk) - (xi + xk) * (beta_Y_squared * xi + xk))
367 / pow(beta_Y_squared * xi + xk, 3);
368 deriv +=
c_Y_ij(i, k, beta, gamma, Y_c_ij) * d2fYikdxidbeta +
dfYikdxi__constxk(x, i, k, beta) * 2 * gamma[i][k] * Y_c_ij[i][k];
382 for (std::size_t k = 0; k < i; k++) {
385 for (std::size_t k = i + 1; k <
N; k++) {
395 for (std::size_t k = 0; k < i; k++) {
398 for (std::size_t k = i + 1; k <
N - 1; k++) {
401 double beta_Y_iN = beta[i][
N - 1], xN = x[
N - 1];
402 deriv += 2 * beta[i][
N - 1] * Y_c_ij[i][
N - 1]
403 * (xN * (x[i] + xN) / (pow(beta_Y_iN, 2) * x[i] + xN)
404 + (1 - beta_Y_iN * beta_Y_iN) * x[i] * xN * xN /
POW2(beta_Y_iN * beta_Y_iN * x[i] + xN));
405 for (std::size_t k = 0; k <
N - 1; ++k) {
406 double beta_Y_kN = beta[k][
N - 1], xk = x[k], beta_Y_kN_squared = beta_Y_kN * beta_Y_kN;
407 deriv += 2 * beta[k][
N - 1] * Y_c_ij[k][
N - 1]
408 * (-xk * (xk + xN) / (beta_Y_kN_squared * xk + xN) + (1 - beta_Y_kN_squared) * xN * xk * xk /
POW2(beta_Y_kN_squared * xk + xN));
421 for (std::size_t k = 0; k < i; k++) {
424 for (std::size_t k = i + 1; k <
N; k++) {
434 for (std::size_t k = 0; k < i; k++) {
437 for (std::size_t k = i + 1; k <
N - 1; k++) {
440 double beta_Y_iN = beta[i][
N - 1], xN = x[
N - 1];
441 d2Yr_dxi2 += 2 *
c_Y_ij(i,
N - 1, beta, gamma, Y_c_ij)
442 * (-(x[i] + xN) / (pow(beta_Y_iN, 2) * x[i] + xN)
443 + (1 - beta_Y_iN * beta_Y_iN)
444 * (xN * xN / pow(beta_Y_iN * beta_Y_iN * x[i] + xN, 2)
445 + ((1 - beta_Y_iN * beta_Y_iN) * x[i] * xN * xN - beta_Y_iN * beta_Y_iN * x[i] * x[i] * xN)
446 / pow(beta_Y_iN * beta_Y_iN * x[i] + xN, 3)));
447 for (std::size_t k = 0; k <
N - 1; ++k) {
448 double beta_Y_kN = beta[k][
N - 1], xk = x[k], beta_Y_kN_squared = beta_Y_kN * beta_Y_kN;
449 d2Yr_dxi2 += 2 *
c_Y_ij(k,
N - 1, beta, gamma, Y_c_ij) * xk * xk * (1 - beta_Y_kN_squared) / pow(beta_Y_kN_squared * xk + xN, 2)
450 * (xN / (beta_Y_kN_squared * xk + xN) - 1);
470 if (j ==
N - 1 || i ==
N - 1) {
479 for (std::size_t k = 0; k <
N - 1; k++) {
494 if (i != j && j != k && k != i) {
496 }
else if (k == i && i != j) {
498 }
else if (k == j && i != j) {
500 }
else if (i == j && i != k) {
504 for (std::size_t m = 0; m < i; m++) {
507 for (std::size_t m = i + 1; m <
N; m++) {
515 for (std::size_t m = 0; m <
N - 1; m++) {
518 if (i != j && j != k && k != i) {
522 }
else if (k == i && i != j) {
525 summer +=
c_Y_ij(i,
N - 1, beta, gamma, Y_c_ij) * (2 *
d3fYijdxidxj2(x, i,
N - 1, beta) -
d3fYijdxi2dxj(x, i,
N - 1, beta));
526 }
else if (k == j && i != j) {
529 summer +=
c_Y_ij(j,
N - 1, beta, gamma, Y_c_ij) * (2 *
d3fYijdxidxj2(x, j,
N - 1, beta) -
d3fYijdxi2dxj(x, j,
N - 1, beta));
530 }
else if (i == j && i != k) {
533 summer +=
c_Y_ij(i,
N - 1, beta, gamma, Y_c_ij) * (2 *
d3fYijdxidxj2(x, i,
N - 1, beta) -
d3fYijdxi2dxj(x, i,
N - 1, beta));
535 for (std::size_t m = 0; m < i; m++) {
538 for (std::size_t m = i + 1; m <
N - 1; m++) {
541 summer +=
c_Y_ij(i,
N - 1, beta, gamma, Y_c_ij)
552 double xk = x[k], xi = x[i], beta_Y = beta[k][i], beta_Y_squared = beta_Y * beta_Y;
553 return xk * (xk + xi) / (beta_Y_squared * xk + xi) + xk * xi / (beta_Y_squared * xk + xi) * (1 - (xk + xi) / (beta_Y_squared * xk + xi));
557 double xk = x[k], xi = x[i], beta_Y = beta[i][k], beta_Y_squared = beta_Y * beta_Y;
558 return xk * (xi + xk) / (beta_Y_squared * xi + xk)
559 + xi * xk / (beta_Y_squared * xi + xk) * (1 - beta_Y_squared * (xi + xk) / (beta_Y_squared * xi + xk));
563 return 2 * beta[i][j] * gamma[i][j] * Y_c[i][j];
566 double xi = x[i], xj = x[j], beta_Y = beta[i][j];
567 return xi * xj * (xi + xj) / (beta_Y * beta_Y * xi + xj);
571 double xi = x[i], xk = x[k], beta_Y = beta[i][k];
572 return 1 / (beta_Y * beta_Y * xi + xk) * (1 - beta_Y * beta_Y * (xi + xk) / (beta_Y * beta_Y * xi + xk))
573 * (2 * xk - xi * xk * 2 * beta_Y * beta_Y / (beta_Y * beta_Y * xi + xk));
577 double xi = x[i], xk = x[k], beta_Y = beta[k][i];
578 return 1 / (beta_Y * beta_Y * xk + xi) * (1 - (xk + xi) / (beta_Y * beta_Y * xk + xi)) * (2 * xk - xk * xi * 2 / (beta_Y * beta_Y * xk + xi));
581 double xi = x[i], xj = x[j], beta_Y = beta[i][j], beta_Y2 = beta_Y * beta_Y;
582 return (xi + xj) / (beta_Y2 * xi + xj) + xj / (beta_Y2 * xi + xj) * (1 - (xi + xj) / (beta_Y2 * xi + xj))
583 + xi / (beta_Y2 * xi + xj) * (1 - beta_Y2 * (xi + xj) / (beta_Y2 * xi + xj))
584 - xi * xj / pow(beta_Y2 * xi + xj, 2) * (1 + beta_Y2 - 2 * beta_Y2 * (xi + xj) / (beta_Y2 * xi + xj));
587 double x_i = x[i], x_j = x[j], beta_Y = beta[i][j], beta_Y2 = beta_Y * beta_Y;
588 double den = pow(beta_Y, 8) * pow(x_i, 4) + 4 * pow(beta_Y, 6) * pow(x_i, 3) * x_j + 6 * pow(beta_Y, 4) * pow(x_i * x_j, 2)
589 + 4 * beta_Y2 * x_i * pow(x_j, 3) + pow(x_j, 4);
590 return -6 * beta_Y2 * x_i * x_j * x_j * (beta_Y2 - 1) / den;
593 double x_i = x[i], x_j = x[j], beta_Y = beta[i][j], beta_Y2 = beta_Y * beta_Y;
594 double den = pow(beta_Y, 8) * pow(x_i, 4) + 4 * pow(beta_Y, 6) * pow(x_i, 3) * x_j + 6 * pow(beta_Y, 4) * pow(x_i * x_j, 2)
595 + 4 * beta_Y2 * x_i * pow(x_j, 3) + pow(x_j, 4);
596 return 6 * beta_Y2 * x_i * x_i * x_j * (beta_Y2 - 1) / den;
600 double x_i = x[i], x_k = x[k], beta_Y = beta[i][k], beta_Y2 = beta_Y * beta_Y;
601 double den = pow(beta_Y, 8) * pow(x_i, 4) + 4 * pow(beta_Y, 6) * pow(x_i, 3) * x_k + 6 * pow(beta_Y, 4) * pow(x_i * x_k, 2)
602 + 4 * beta_Y2 * x_i * pow(x_k, 3) + pow(x_k, 4);
603 return 6 * beta_Y2 * x_k * x_k * x_k * (beta_Y2 - 1) / den;
607 double x_i = x[i], x_k = x[k], beta_Y = beta[k][i], beta_Y2 = beta_Y * beta_Y;
608 return 6 * beta_Y2 * x_k * x_k * x_k * (1 - beta_Y2) / pow(beta_Y2 * x_k + x_i, 4);