Skip to content

Commit

Permalink
improvements to kinetic API
Browse files Browse the repository at this point in the history
  • Loading branch information
@yoelcortes
yoelcortes committed Jun 4, 2024
1 parent 5e4ef7d commit 800b77f
Show file tree
Hide file tree
Showing 6 changed files with 275 additions and 241 deletions.
14 changes: 7 additions & 7 deletions tests/test_reaction.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -223,7 +223,7 @@ def test_reactive_phase_equilibrium_no_kinetics():
stream = tmo.Stream(
H2O=1, Ethanol=5, LacticAcid=1
)
stream.vle(T=360, P=101325, liquid_reaction=rxn)
stream.vle(T=360, P=101325, liquid_conversion=rxn)
assert_allclose(
stream.imol['g'],
[0.0,
Expand All @@ -239,7 +239,7 @@ def test_reactive_phase_equilibrium_no_kinetics():
1.8234109156732896]
)

stream.vle(T=340, P=101325, liquid_reaction=rxn)
stream.vle(T=340, P=101325, liquid_conversion=rxn)
assert_allclose(
stream.imol['l'],
stream.mol
Expand All @@ -248,7 +248,7 @@ def test_reactive_phase_equilibrium_no_kinetics():
stream.imol['g'],
0
)
stream.vle(T=450, P=101325, liquid_reaction=rxn)
stream.vle(T=450, P=101325, liquid_conversion=rxn)
assert_allclose(
stream.imol['g'],
stream.mol
Expand Down Expand Up @@ -297,7 +297,7 @@ def rate(self, stream):
)
T = 360
P = 101325
stream.vle(T=T, P=P, liquid_reaction=rxn)
stream.vle(T=T, P=P, liquid_conversion=rxn)
assert_allclose(
stream.imol['l'],
[0.026722998037919946,
Expand All @@ -321,7 +321,7 @@ def rate(self, stream):
stream = tmo.Stream(
H2O=2, Ethanol=5, LacticAcid=1, T=T,
)
stream.vle(V=V, P=P, liquid_reaction=rxn)
stream.vle(V=V, P=P, liquid_conversion=rxn)
assert_allclose(
stream.imol['l'],
[0.026426291229780553,
Expand All @@ -343,7 +343,7 @@ def rate(self, stream):
stream = tmo.Stream(
H2O=2, Ethanol=5, LacticAcid=1, T=T,
)
stream.vle(V=V, T=T, liquid_reaction=rxn)
stream.vle(V=V, T=T, liquid_conversion=rxn)
assert_allclose(
stream.imol['l'],
[0.026556271540719167,
Expand All @@ -362,7 +362,7 @@ def rate(self, stream):
stream = tmo.Stream(
H2O=2, Ethanol=5, LacticAcid=1, T=T,
)
stream.vle(H=H, P=P, liquid_reaction=rxn)
stream.vle(H=H, P=P, liquid_conversion=rxn)
assert_allclose(
stream.imol['l'],
[0.026722998039327987,
Expand Down
32 changes: 31 additions & 1 deletion thermosteam/_stream.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -60,12 +60,39 @@ def __init__(self, stream, original, temporary):
self.temporary = temporary

def __enter__(self):
self.stream._phase._phase = self.temporary
stream = self.stream
stream._phase._phase = self.temporary
return stream

def __exit__(self, type, exception, traceback):
self.stream._phase._phase = self.original
if exception: raise exception


class TemporaryStream:
__slots__ = ('stream', 'data', 'flow', 'T', 'P', 'phase')

def __init__(self, stream, flow, T, P, phase):
self.stream = stream
self.data = stream.get_data()
self.flow = flow
self.T = T
self.P = P
self.phase = phase

def __enter__(self):
stream = self.stream
self.data = stream.get_data()
if self.flow is not None: stream.imol.data[:] = self.flow
if self.T is not None: stream.T = self.T
if self.P is not None: stream.P = self.P
return stream

def __exit__(self, type, exception, traceback):
self.stream.set_data(self.data)
if exception: raise exception


class Equations:
__slots__ = ('material', 'energy')

Expand Down Expand Up @@ -361,6 +388,9 @@ def __init__(self, ID: Optional[str]='',
data = self._imol.data
self.phases = [j for i, j in enumerate(self.phases) if data[i].any()]

def temporary(self, flow=None, T=None, P=None, phase=None):
return TemporaryStream(self, flow, T, P, phase)

def temporary_phase(self, phase):
return TemporaryPhase(self, self.phase, phase)

Expand Down
32 changes: 16 additions & 16 deletions thermosteam/equilibrium/bubble_point.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -144,9 +144,9 @@ def _P_error(self, P, T, z_Psat_gamma, Psats, y):
y[:] = solve_y(y_phi, self.phi, T, P, y)
return 1. - y.sum()

def _T_error_reactive(self, T, P, z, dz, y, x, liquid_reaction):
def _T_error_reactive(self, T, P, z, dz, y, x, liquid_conversion):
if T <= 0: raise InfeasibleRegion('negative temperature')
dz[:] = liquid_reaction.conversion(z, T, P, 'l')
dz[:] = liquid_conversion(z, T, P, 'l')
x[:] = z + dz
x /= x.sum()
Psats = np.array([i(T) for i in self.Psats], dtype=float)
Expand All @@ -157,9 +157,9 @@ def _T_error_reactive(self, T, P, z, dz, y, x, liquid_reaction):
y[:] = solve_y(y_phi, self.phi, T, P, y)
return 1. - y.sum()

def _P_error_reactive(self, P, T, Psats, z, dz, y, x, liquid_reaction):
def _P_error_reactive(self, P, T, Psats, z, dz, y, x, liquid_conversion):
if P <= 0: raise InfeasibleRegion('negative pressure')
dz[:] = liquid_reaction.conversion(z, T, P, 'l')
dz[:] = liquid_conversion(z, T, P, 'l')
x[:] = z + dz
x /= x.sum()
z_Psat_gamma = x * Psats * self.gamma(x, T)
Expand Down Expand Up @@ -192,21 +192,21 @@ def _Py_ideal(self, z_Psat_gamma_pcf):
y = z_Psat_gamma_pcf / P
return P, y

def __call__(self, z, *, T=None, P=None, liquid_reaction=None):
def __call__(self, z, *, T=None, P=None, liquid_conversion=None):
z = np.asarray(z, float)
if T:
if P: raise ValueError("may specify either T or P, not both")
P, *args = self.solve_Py(z, T, liquid_reaction)
P, *args = self.solve_Py(z, T, liquid_conversion)
elif P:
T, *args = self.solve_Ty(z, P, liquid_reaction)
T, *args = self.solve_Ty(z, P, liquid_conversion)
else:
raise ValueError("must specify either T or P")
if liquid_reaction:
if liquid_conversion:
return ReactiveBubblePointValues(T, P, self.IDs, z, *args)
else:
return BubblePointValues(T, P, self.IDs, z, *args)

def solve_Ty(self, z, P, liquid_reaction=None):
def solve_Ty(self, z, P, liquid_conversion=None):
"""
Bubble point at given composition and pressure.
Expand Down Expand Up @@ -244,7 +244,7 @@ def solve_Ty(self, z, P, liquid_reaction=None):
T = chemical.Tsat(P, check_validity=False) if P <= chemical.Pc else chemical.Tc
y = z.copy()
return T, fn.normalize(y)
elif liquid_reaction is None:
elif liquid_conversion is None:
f = self._T_error
z_norm = z / z.sum()
z_over_P = z/P
Expand All @@ -268,7 +268,7 @@ def solve_Ty(self, z, P, liquid_reaction=None):
dz = z_norm.copy()
z_over_P = z / P
T_guess, y = self._Ty_ideal(z_over_P)
args = (P, z_norm, dz, y, x, liquid_reaction)
args = (P, z_norm, dz, y, x, liquid_conversion)
try:
T = flx.aitken_secant(f, T_guess, T_guess + 1e-3,
self.T_tol, 5e-12, args,
Expand All @@ -281,7 +281,7 @@ def solve_Ty(self, z, P, liquid_reaction=None):
checkiter=False, checkbounds=False)
return T, dz, fn.normalize(y), x

def solve_Py(self, z, T, liquid_reaction=None):
def solve_Py(self, z, T, liquid_conversion=None):
"""
Bubble point at given composition and temperature.
Expand Down Expand Up @@ -319,7 +319,7 @@ def solve_Py(self, z, T, liquid_reaction=None):
P = chemical.Psat(T) if T <= chemical.Tc else chemical.Pc
y = z.copy()
return P, fn.normalize(y)
elif liquid_reaction is None:
elif liquid_conversion is None:
if T > self.Tmax: T = self.Tmax
elif T < self.Tmin: T = self.Tmin
Psats = np.array([i(T) for i in self.Psats])
Expand All @@ -346,7 +346,7 @@ def solve_Py(self, z, T, liquid_reaction=None):
dz = z_norm.copy()
z_Psat_gamma = z * Psats * self.gamma(z_norm, T)
P_guess, y = self._Py_ideal(z_Psat_gamma)
args = (T, Psats, z_norm, dz, y, x, liquid_reaction)
args = (T, Psats, z_norm, dz, y, x, liquid_conversion)
try:
P = flx.aitken_secant(f, P_guess, P_guess-1, self.P_tol, 1e-9,
args, checkiter=False)
Expand Down Expand Up @@ -400,7 +400,7 @@ def __init__(self, chemicals=(), flasher=None):

__call__ = BubblePoint.__call__

def solve_Ty(self, z, P, liquid_reaction=None):
def solve_Ty(self, z, P, liquid_conversion=None):
"""
Bubble point at given composition and pressure.
Expand Down Expand Up @@ -442,7 +442,7 @@ def solve_Ty(self, z, P, liquid_reaction=None):
T = results.T
return T, fn.normalize(y)

def solve_Py(self, z, T, liquid_reaction=None):
def solve_Py(self, z, T, liquid_conversion=None):
"""
Bubble point at given composition and temperature.
Expand Down
28 changes: 14 additions & 14 deletions thermosteam/equilibrium/dew_point.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -151,9 +151,9 @@ def _T_error(self, T, P, z_norm, zP, x):
x[:] = self._solve_x(x_gamma, T, P, x)
return 1 - x.sum()

def _T_error_reactive(self, T, P, z, dz, y, x, gas_reaction):
def _T_error_reactive(self, T, P, z, dz, y, x, gas_conversion):
if T <= 0: raise InfeasibleRegion('negative temperature')
dz[:] = gas_reaction.conversion(z, T, P, 'g')
dz[:] = gas_conversion(z, T, P, 'g')
y[:] = z + dz
y /= y.sum()
Psats = np.array([i(T) for i in self.Psats])
Expand All @@ -176,9 +176,9 @@ def _P_error(self, P, T, z_norm, z_over_Psats, Psats, x):
x[:] = self._solve_x(x_gamma, T, P, x)
return 1 - x.sum()

def _P_error_reactive(self, P, T, Psats, z, dz, y, x, gas_reaction):
def _P_error_reactive(self, P, T, Psats, z, dz, y, x, gas_conversion):
if P <= 0: raise InfeasibleRegion('negative pressure')
dz[:] = gas_reaction.conversion(z, T, P, 'g')
dz[:] = gas_conversion(z, T, P, 'g')
y[:] = z + dz
y /= y.sum()
x_gamma = y / Psats * P * self.phi(y, T, P) / self.pcf(T, P, Psats)
Expand All @@ -205,21 +205,21 @@ def _Px_ideal(self, z_over_Psats):
x = z_over_Psats * P
return P, x

def __call__(self, z, *, T=None, P=None, gas_reaction=None):
def __call__(self, z, *, T=None, P=None, gas_conversion=None):
z = np.asarray(z, float)
if T:
if P: raise ValueError("may specify either T or P, not both")
P, *args = self.solve_Px(z, T, gas_reaction)
P, *args = self.solve_Px(z, T, gas_conversion)
elif P:
T, *args = self.solve_Tx(z, P, gas_reaction)
T, *args = self.solve_Tx(z, P, gas_conversion)
else:
raise ValueError("must specify either T or P")
if gas_reaction:
if gas_conversion:
return ReactiveDewPointValues(T, P, self.IDs, z, *args)
else:
return DewPointValues(T, P, self.IDs, z, *args)

def solve_Tx(self, z, P, gas_reaction=None):
def solve_Tx(self, z, P, gas_conversion=None):
"""
Dew point given composition and pressure.
Expand Down Expand Up @@ -257,7 +257,7 @@ def solve_Tx(self, z, P, gas_reaction=None):
T = chemical.Tsat(P, check_validity=False) if P <= chemical.Pc else chemical.Tc
x = z.copy()
return T, fn.normalize(x)
elif gas_reaction is None:
elif gas_conversion is None:
f = self._T_error
z_norm = z/z.sum()
zP = z * P
Expand All @@ -282,7 +282,7 @@ def solve_Tx(self, z, P, gas_reaction=None):
dz = z_norm.copy()
zP = z * P
T_guess, y = self._Tx_ideal(zP)
args = (P, z_norm, dz, y, x, gas_reaction)
args = (P, z_norm, dz, y, x, gas_conversion)
try:
T = flx.aitken_secant(f, T_guess, T_guess + 1e-3,
self.T_tol, 5e-12, args,
Expand All @@ -295,7 +295,7 @@ def solve_Tx(self, z, P, gas_reaction=None):
checkiter=False, checkbounds=False)
return T, dz, fn.normalize(y), x

def solve_Px(self, z, T, gas_reaction=None):
def solve_Px(self, z, T, gas_conversion=None):
"""
Dew point given composition and temperature.
Expand Down Expand Up @@ -333,7 +333,7 @@ def solve_Px(self, z, T, gas_reaction=None):
P = chemical.Psat(T) if T <= chemical.Tc else chemical.Pc
x = z.copy()
return P, fn.normalize(x)
elif gas_reaction is None:
elif gas_conversion is None:
z_norm = z / z.sum()
Psats = np.array([i(T) for i in self.Psats], dtype=float)
z_over_Psats = z / Psats
Expand All @@ -359,7 +359,7 @@ def solve_Px(self, z, T, gas_reaction=None):
Psats = np.array([i(T) for i in self.Psats], dtype=float)
z_over_Psats = z / Psats
P_guess, x = self._Px_ideal(z_over_Psats)
args = (T, Psats, z_norm, dz, y, x, gas_reaction)
args = (T, Psats, z_norm, dz, y, x, gas_conversion)
try:
P = flx.aitken_secant(f, P_guess, P_guess-1, self.P_tol, 1e-9,
args, checkiter=False)
Expand Down
Loading

0 comments on commit 800b77f

Please sign in to comment.