Skip to content

Commit d21633a

Browse files
committed
Improvements to NRTL
1 parent 3c5bc40 commit d21633a

File tree

2 files changed

+71
-178
lines changed

2 files changed

+71
-178
lines changed

tests/test_nrtl.py

Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -32,7 +32,7 @@
3232
import numpy as np
3333
from fluids.numerics import jacobian, hessian, derivative, normalize, assert_close, assert_close1d, assert_close2d
3434

35-
def test_NRTL():
35+
def test_NRTL_gammas():
3636
# P05.01b VLE Behavior of Ethanol - Water Using NRTL
3737
gammas = NRTL_gammas([0.252, 0.748], [[0, -0.178], [1.963, 0]], [[0, 0.2974],[.2974, 0]])
3838
assert_close1d(gammas, [1.9363183763514304, 1.1537609663170014])
@@ -96,7 +96,7 @@ def NRTL2(xs, taus, alpha):
9696
gammas_expect = [2.503204848288857, 2.910723989902569, 2.2547951278295497, 2.9933258413917154, 2.694165187439594]
9797
assert_close1d(gammas, gammas_expect)
9898

99-
def test_NRTL_10():
99+
def test_NRTL_gammas_10():
100100
# ten component
101101
# m = Mixture(['water', 'ethanol', 'methanol', '1-pentanol', '2-pentanol', '3-pentanol',
102102
# '1-decanol', '2-decanol', '3-decanol', '4-decanol'],
@@ -129,15 +129,15 @@ def test_NRTL_10():
129129
gammas_expect = [1.1600804309840225, 1.0892286716705042, 1.0384940848807305, 0.9836770920034531, 0.9836770920034531, 0.9836770920034531, 0.9836770920034531, 0.9836770920034531, 0.9836770920034531, 0.9836770920034531]
130130
assert_close1d(gammas, gammas_expect)
131131

132-
def NRTL_44():
132+
def NRTL_gammas_44():
133133
N = 44
134134
taus = [[random() for i in range(N)] for j in range(N)]
135135
alphas = [[random() for i in range(N)] for j in range(N)]
136136
xs = normalize([random() for i in range(N)])
137137
gammas = NRTL_gammas(xs=xs, taus=taus, alphas=alphas)
138138

139139

140-
def NRTL_200():
140+
def NRTL_gammas_200():
141141
# ten component
142142
# Takes 40 ms - not a great idea
143143
N = 200
@@ -175,8 +175,8 @@ def make_taus(N):
175175
return data
176176

177177

178-
def test_madeup_20():
179-
N = 10
178+
def test_madeup_NRTL():
179+
N = 6
180180
alphas = make_alphas(N)
181181
taus = make_taus(N)
182182
xs = normalize([random() for i in range(N)])

thermo/nrtl.py

Lines changed: 65 additions & 172 deletions
Original file line numberDiff line numberDiff line change
@@ -692,22 +692,12 @@ def gammas(self):
692692

693693
xj_Gs_jis_inv, xj_Gs_taus_jis = self.xj_Gs_jis_inv(), self.xj_Gs_taus_jis()
694694

695-
self._gammas = gammas = nrtl_gammas(xs, N, Gs, taus, xj_Gs_jis_inv, xj_Gs_taus_jis)
696-
697-
# t0s = [xs[j]*xj_Gs_jis_inv[j] for j in range(N)]
698-
# t1s = [xj_Gs_taus_jis[j]*xj_Gs_jis_inv[j] for j in range(N)]
699-
#
700-
# for i in range(N):
701-
# tot = xj_Gs_taus_jis[i]*xj_Gs_jis_inv[i]
702-
# Gsi = Gs[i]
703-
# tausi = taus[i]
704-
# for j in range(N):
705-
# # Possible to factor out some terms which depend on j only; or to index taus, Gs separately
706-
# tot += t0s[j]*Gsi[j]*(tausi[j] - t1s[j])
707-
#
708-
# gammas.append(exp(tot))
695+
if self.scalar:
696+
gammas = [0.0]*self.N
697+
else:
698+
gammas = zeros(self.N)
709699

710-
# self._gammas = gammas = NRTL_gammas(xs=self.xs, taus=taus, alphas=alphas)
700+
self._gammas = nrtl_gammas(xs, N, Gs, taus, xj_Gs_jis_inv, xj_Gs_taus_jis, gammas)
711701
return gammas
712702

713703

@@ -738,35 +728,13 @@ def taus(self):
738728
F = self.tau_coeffs_F
739729
G = self.tau_coeffs_G
740730
H = self.tau_coeffs_H
741-
742-
# tau_coeffs_A = self.tau_coeffs_A
743-
# tau_coeffs_B = self.tau_coeffs_B
744-
# tau_coeffs_E = self.tau_coeffs_E
745-
# tau_coeffs_F = self.tau_coeffs_F
746-
# tau_coeffs_G = self.tau_coeffs_G
747-
# tau_coeffs_H = self.tau_coeffs_H
748-
749-
750731
T, N = self.T, self.N
751-
# T2 = T*T
752-
# Tinv = 1.0/T
753-
# T2inv = Tinv*Tinv
754-
# logT = log(T)
755-
#
756-
# initialize the matrix to be A
757-
self._taus = taus = nrtl_taus(T, N, A, B, E, F, G, H)
758-
# self._taus = taus = [list(l) for l in tau_coeffs_A]
759-
# for i in range(N):
760-
# tau_coeffs_Bi = tau_coeffs_B[i]
761-
# tau_coeffs_Ei = tau_coeffs_E[i]
762-
# tau_coeffs_Fi = tau_coeffs_F[i]
763-
# tau_coeffs_Gi = tau_coeffs_G[i]
764-
# tau_coeffs_Hi = tau_coeffs_H[i]
765-
# tausi = taus[i]
766-
# for j in range(N):
767-
# tausi[j] += (tau_coeffs_Bi[j]*Tinv + tau_coeffs_Ei[j]*logT
768-
# + tau_coeffs_Fi[j]*T + tau_coeffs_Gi[j]*T2inv
769-
# + tau_coeffs_Hi[j]*T2)
732+
if self.scalar:
733+
taus = [[0.0]*N for _ in range(N)]
734+
else:
735+
taus = zeros((N, N))
736+
737+
self._taus = nrtl_taus(T, N, A, B, E, F, G, H, taus)
770738
return taus
771739

772740
def dtaus_dT(self):
@@ -797,24 +765,11 @@ def dtaus_dT(self):
797765
G = self.tau_coeffs_G
798766
H = self.tau_coeffs_H
799767
T, N = self.T, self.N
800-
801-
# Tinv = 1.0/T
802-
# nT2inv = -Tinv*Tinv
803-
# n2T3inv = 2.0*nT2inv*Tinv
804-
# T2 = T + T
805-
806-
# self._dtaus_dT = dtaus_dT = [list(l) for l in tau_coeffs_F]
807-
self._dtaus_dT = dtaus_dT = nrtl_dtaus_dT(T, N, B, E, F, G, H)
808-
# for i in range(N):
809-
# tau_coeffs_Bi = tau_coeffs_B[i]
810-
# tau_coeffs_Ei = tau_coeffs_E[i]
811-
# tau_coeffs_Fi = tau_coeffs_F[i]
812-
# tau_coeffs_Gi = tau_coeffs_G[i]
813-
# tau_coeffs_Hi = tau_coeffs_H[i]
814-
# dtaus_dTi = dtaus_dT[i]
815-
# for j in range(N):
816-
# dtaus_dTi[j] += (nT2inv*tau_coeffs_Bi[j] + Tinv*tau_coeffs_Ei[j]
817-
# + n2T3inv*tau_coeffs_Gi[j] + T2*tau_coeffs_Hi[j])
768+
if self.scalar:
769+
dtaus_dT = [[0.0]*N for _ in range(N)]
770+
else:
771+
dtaus_dT = zeros((N, N))
772+
self._dtaus_dT = nrtl_dtaus_dT(T, N, B, E, F, G, H, dtaus_dT)
818773

819774
return dtaus_dT
820775

@@ -845,24 +800,12 @@ def d2taus_dT2(self):
845800
H = self.tau_coeffs_H
846801
T, N = self.T, self.N
847802

848-
# self._d2taus_dT2 = d2taus_dT2 = [[h + h for h in l] for l in tau_coeffs_H]
849-
self._d2taus_dT2 = d2taus_dT2 = nrtl_d2taus_dT2(T, N, B, E, G, H)
803+
if self.scalar:
804+
d2taus_dT2 = [[0.0]*N for _ in range(N)]
805+
else:
806+
d2taus_dT2 = zeros((N, N))
850807

851-
# Tinv = 1.0/T
852-
# Tinv2 = Tinv*Tinv
853-
#
854-
# T3inv2 = 2.0*(Tinv2*Tinv)
855-
# nT2inv = -Tinv*Tinv
856-
# T4inv6 = 6.0*(Tinv2*Tinv2)
857-
# for i in range(N):
858-
# tau_coeffs_Bi = tau_coeffs_B[i]
859-
# tau_coeffs_Ei = tau_coeffs_E[i]
860-
# tau_coeffs_Gi = tau_coeffs_G[i]
861-
# d2taus_dT2i = d2taus_dT2[i]
862-
# for j in range(N):
863-
# d2taus_dT2i[j] += (T3inv2*tau_coeffs_Bi[j]
864-
# + nT2inv*tau_coeffs_Ei[j]
865-
# + T4inv6*tau_coeffs_Gi[j])
808+
self._d2taus_dT2 = nrtl_d2taus_dT2(T, N, B, E, G, H, d2taus_dT2)
866809
return d2taus_dT2
867810

868811
def d3taus_dT3(self):
@@ -890,25 +833,11 @@ def d3taus_dT3(self):
890833
E = self.tau_coeffs_E
891834
G = self.tau_coeffs_G
892835
T, N = self.T, self.N
893-
# self._d3taus_dT3 = d3taus_dT3 = [[0.0]*N for i in range(N)]
894-
self._d3taus_dT3 = d3taus_dT3 = nrtl_d3taus_dT3(T, N, B, E, G)
895-
896-
# Tinv = 1.0/T
897-
# T2inv = Tinv*Tinv
898-
#
899-
# nT4inv6 = -6.0*T2inv*T2inv
900-
# T3inv2 = 2.0*T2inv*Tinv
901-
# T5inv24 = -24.0*(T2inv*T2inv*Tinv)
902-
903-
# for i in range(N):
904-
# tau_coeffs_Bi = tau_coeffs_B[i]
905-
# tau_coeffs_Ei = tau_coeffs_E[i]
906-
# tau_coeffs_Gi = tau_coeffs_G[i]
907-
# d3taus_dT3i = d3taus_dT3[i]
908-
# for j in range(N):
909-
# d3taus_dT3i[j] = (nT4inv6*tau_coeffs_Bi[j]
910-
# + T3inv2*tau_coeffs_Ei[j]
911-
# + T5inv24*tau_coeffs_Gi[j])
836+
if self.scalar:
837+
d3taus_dT3 = [[0.0]*N for _ in range(N)]
838+
else:
839+
d3taus_dT3 = zeros((N, N))
840+
self._d3taus_dT3 = nrtl_d3taus_dT3(T, N, B, E, G, d3taus_dT3)
912841
return d3taus_dT3
913842

914843
def alphas(self):
@@ -944,15 +873,14 @@ def alphas(self):
944873
return self._alphas
945874
except AttributeError:
946875
pass
947-
# T, N = self.T, self.N
948-
# alpha_coeffs_c, alpha_coeffs_d = self.alpha_coeffs_c, self.alpha_coeffs_d
876+
N = self.N
949877

950-
self._alphas = alphas = nrtl_alphas(self.T, self.N, self.alpha_coeffs_c, self.alpha_coeffs_d)
951-
# for i in range(N):
952-
# alpha_coeffs_ci = alpha_coeffs_c[i]
953-
# alpha_coeffs_di = alpha_coeffs_d[i]
954-
# alphas.append([alpha_coeffs_ci[j] + alpha_coeffs_di[j]*T for j in range(N)])
878+
if self.scalar:
879+
alphas = [[0.0]*N for _ in range(N)]
880+
else:
881+
alphas = zeros((N, N))
955882

883+
self._alphas = nrtl_alphas(self.T, N, self.alpha_coeffs_c, self.alpha_coeffs_d, alphas)
956884
return alphas
957885

958886
def dalphas_dT(self):
@@ -990,11 +918,12 @@ def Gs(self):
990918
taus = self.taus()
991919
N = self.N
992920

993-
self._Gs = Gs = nrtl_Gs(N, alphas, taus)
994-
# for i in range(N):
995-
# alphasi = alphas[i]
996-
# tausi = taus[i]
997-
# Gs.append([exp(-alphasi[j]*tausi[j]) for j in range(N)])
921+
if self.scalar:
922+
Gs = [[0.0]*N for _ in range(N)]
923+
else:
924+
Gs = zeros((N, N))
925+
926+
self._Gs = nrtl_Gs(N, alphas, taus, Gs)
998927
return Gs
999928

1000929
def dGs_dT(self):
@@ -1031,16 +960,12 @@ def dGs_dT(self):
1031960
Gs = self.Gs()
1032961
N = self.N
1033962

1034-
self._dGs_dT = dGs_dT = nrtl_dGs_dT(N, alphas, dalphas_dT, taus, dtaus_dT, Gs)
1035-
# for i in range(N):
1036-
# alphasi = alphas[i]
1037-
# tausi = taus[i]
1038-
# dalphasi = dalphas_dT[i]
1039-
# dtausi = dtaus_dT[i]
1040-
# Gsi = Gs[i]
1041-
#
1042-
# dGs_dT.append([(-alphasi[j]*dtausi[j] - tausi[j]*dalphasi[j])*Gsi[j]
1043-
# for j in range(N)])
963+
if self.scalar:
964+
dGs_dT = [[0.0]*N for _ in range(N)]
965+
else:
966+
dGs_dT = zeros((N, N))
967+
968+
self._dGs_dT = nrtl_dGs_dT(N, alphas, dalphas_dT, taus, dtaus_dT, Gs, dGs_dT)
1044969
return dGs_dT
1045970

1046971
def d2Gs_dT2(self):
@@ -1083,23 +1008,12 @@ def d2Gs_dT2(self):
10831008
Gs = self.Gs()
10841009
N = self.N
10851010

1086-
self._d2Gs_dT2 = d2Gs_dT2 = nrtl_d2Gs_dT2(N, alphas, dalphas_dT, taus, dtaus_dT, d2taus_dT2, Gs)
1087-
# nrtl_d2Gs_dT2(N, alphas, dalphas_dT, taus, dtaus_dT, d2taus_dT2, Gs
1088-
# []
1089-
# for i in range(N):
1090-
# alphasi = alphas[i]
1091-
# tausi = taus[i]
1092-
# dalphasi = dalphas_dT[i]
1093-
# dtausi = dtaus_dT[i]
1094-
# d2taus_dT2i = d2taus_dT2[i]
1095-
# Gsi = Gs[i]
1096-
#
1097-
# d2Gs_dT2_row = []
1098-
# for j in range(N):
1099-
# t1 = alphasi[j]*dtausi[j] + tausi[j]*dalphasi[j]
1100-
# d2Gs_dT2_row.append((t1*t1 - alphasi[j]*d2taus_dT2i[j]
1101-
# - 2.0*dalphasi[j]*dtausi[j])*Gsi[j])
1102-
# d2Gs_dT2.append(d2Gs_dT2_row)
1011+
if self.scalar:
1012+
d2Gs_dT2 = [[0.0]*N for _ in range(N)]
1013+
else:
1014+
d2Gs_dT2 = zeros((N, N))
1015+
1016+
self._d2Gs_dT2 = nrtl_d2Gs_dT2(N, alphas, dalphas_dT, taus, dtaus_dT, d2taus_dT2, Gs, d2Gs_dT2)
11031017
return d2Gs_dT2
11041018

11051019
def d3Gs_dT3(self):
@@ -1146,27 +1060,12 @@ def d3Gs_dT3(self):
11461060
d3taus_dT3 = self.d3taus_dT3()
11471061
Gs = self.Gs()
11481062

1149-
self._d3Gs_dT3 = d3Gs_dT3 = nrtl_d3Gs_dT3(N, alphas, dalphas_dT, taus, dtaus_dT, d2taus_dT2, d3taus_dT3, Gs)
1150-
# for i in range(N):
1151-
# alphasi = alphas[i]
1152-
# tausi = taus[i]
1153-
# dalphasi = dalphas_dT[i]
1154-
# dtaus_dTi = dtaus_dT[i]
1155-
# d2taus_dT2i = d2taus_dT2[i]
1156-
# d3taus_dT3i = d3taus_dT3[i]
1157-
# Gsi = Gs[i]
1158-
# d3Gs_dT3_row = []
1159-
# for j in range(N):
1160-
# x0 = alphasi[j]
1161-
# x1 = tausi[j]
1162-
# x2 = dalphasi[j]
1163-
#
1164-
# x3 = d2taus_dT2i[j]
1165-
# x4 = dtaus_dTi[j]
1166-
# x5 = x0*x4 + x1*x2
1167-
# v = Gsi[j]*(-x0*d3taus_dT3i[j] - 3.0*x2*x3 - x5*x5*x5 + 3.0*x5*(x0*x3 + 2.0*x2*x4))
1168-
# d3Gs_dT3_row.append(v)
1169-
# d3Gs_dT3.append(d3Gs_dT3_row)
1063+
if self.scalar:
1064+
d3Gs_dT3 = [[0.0]*N for _ in range(N)]
1065+
else:
1066+
d3Gs_dT3 = zeros((N, N))
1067+
1068+
self._d3Gs_dT3 = nrtl_d3Gs_dT3(N, alphas, dalphas_dT, taus, dtaus_dT, d2taus_dT2, d3taus_dT3, Gs, d3Gs_dT3)
11701069
return d3Gs_dT3
11711070

11721071

@@ -1187,23 +1086,17 @@ def xj_Gs_jis(self):
11871086
taus = self.taus()
11881087

11891088
xs, N = self.xs, self.N
1190-
_xj_Gs_jis, _xj_Gs_taus_jis = nrtl_xj_Gs_jis_and_Gs_taus_jis(N, xs, Gs, taus)
1191-
if not self.scalar and type(_xj_Gs_jis) is list:
1192-
_xj_Gs_jis, _xj_Gs_taus_jis = array(_xj_Gs_jis), array(_xj_Gs_taus_jis)
1089+
if self.scalar:
1090+
_xj_Gs_jis = [0.0]*N
1091+
_xj_Gs_taus_jis = [0.0]*N
1092+
else:
1093+
_xj_Gs_jis = zeros(N)
1094+
_xj_Gs_taus_jis = zeros(N)
1095+
1096+
nrtl_xj_Gs_jis_and_Gs_taus_jis(N, xs, Gs, taus, _xj_Gs_jis, _xj_Gs_taus_jis)
11931097

11941098
self._xj_Gs_jis, self._xj_Gs_taus_jis = _xj_Gs_jis, _xj_Gs_taus_jis
1195-
# self._xj_Gs_jis = xj_Gs_jis = []
1196-
# self._xj_Gs_taus_jis = xj_Gs_taus_jis = []
1197-
# for i in range(N):
1198-
# tot1 = 0.0
1199-
# tot2 = 0.0
1200-
# for j in range(N):
1201-
# xjGji = xs[j]*Gs[j][i]
1202-
# tot1 += xjGji#xs[j]*Gs[j][i]
1203-
# tot2 += xjGji*taus[j][i]
1204-
# xj_Gs_jis.append(tot1)
1205-
# xj_Gs_taus_jis.append(tot2)
1206-
return self._xj_Gs_jis
1099+
return _xj_Gs_jis
12071100

12081101
def xj_Gs_jis_inv(self):
12091102
try:

0 commit comments

Comments
 (0)