Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
34 changes: 20 additions & 14 deletions src/autochem/rate/data.py
Original file line number Diff line number Diff line change
Expand Up @@ -648,7 +648,10 @@ def fit(
if not validate or bad_fit == "fill":
try:
with np.errstate(all="raise"):
obj(T=T)
vals = obj(T=T)
if not np.all(np.isfinite(vals)):
msg = "Fitted rate gives non-finite values over the input temperature range."
raise FloatingPointError(msg)
except FloatingPointError as e:
if A_fill is not None:
return cls(order=order, A=A_fill, b=0, E=0)
Expand Down Expand Up @@ -917,22 +920,25 @@ def fit( # noqa: PLR0913
(Options: "fill", replace with the fill value, or numpy.seterr options)
:return: Rate fit
"""
bad_fits = [
"fill" if np.any(np.isclose(p, bad_fit_fill_pressures)) else bad_fit
for p in Ps
]
k_data_fits = [
ArrheniusRateFit.fit(
k_fits = []
for k_data, P in zip(np.transpose(k_data), Ps, strict=True):
if validate:
print(f"Fitting Arrhenius rate for {P = }")

bad_fit = (
"fill" if np.any(np.isclose(P, bad_fit_fill_pressures)) else bad_fit
)
k_fit = ArrheniusRateFit.fit(
Ts=Ts,
ks=ks,
ks=k_data,
A_fill=A_fill,
bad_fit=f, # pyright: ignore[reportArgumentType]
bad_fit=bad_fit,
validate=validate,
order=order,
units=units,
)
for ks, f in zip(np.transpose(k_data), bad_fits, strict=True)
]
k_fits.append(k_fit)

k_high_fit = None
if k_high is not None:
k_high_fit = ArrheniusRateFit.fit(
Expand All @@ -943,9 +949,9 @@ def fit( # noqa: PLR0913

return cls(
order=order,
As=[f.A for f in k_data_fits],
bs=[f.b for f in k_data_fits],
Es=[f.E for f in k_data_fits],
As=[f.A for f in k_fits],
bs=[f.b for f in k_fits],
Es=[f.E for f in k_fits],
Ps=Ps, # pyright: ignore[reportArgumentType]
)

Expand Down
8 changes: 3 additions & 5 deletions src/autochem/tests/test__therm.py
Original file line number Diff line number Diff line change
Expand Up @@ -81,10 +81,8 @@ def test__from_chemkin_string(name, spc_str0):
# Plot
therm.display(spc)
therm.display(
spc,
label="original",
others=[spc_times_2, spc_divided_by_2],
others_labels=["doubled", "halved"],
[spc, spc_times_2, spc_divided_by_2],
label=["original","doubled", "halved"]
)


Expand All @@ -105,7 +103,7 @@ def test__fit(name, spc_str0):
units={"energy": "kcal"},
)
spc_fit = therm.fit(spc)
therm.display(spc, label="data", others=[spc_fit], others_labels=["fit"])
therm.display([spc, spc_fit], label=["data", "fit"])


if __name__ == "__main__":
Expand Down
35 changes: 26 additions & 9 deletions src/autochem/therm/_species.py
Original file line number Diff line number Diff line change
@@ -1,8 +1,9 @@
"""Thermodynamic data."""

import datetime
import re
from collections.abc import Sequence
from typing import ClassVar, Literal
from typing import ClassVar, Literal, Self

import altair as alt
import numpy as np
Expand All @@ -23,6 +24,17 @@ class Species(Scalable):
name: str
therm: Therm_

def racemize(
self, enant_suffixes: tuple[str, str] = ("0", "1"), rac_suffix: str = "R"
) -> Self:
"""Racemize the thermodynamic functions for a chiral species."""
suffix0, suffix1 = enant_suffixes
pattern = re.compile(rf"(?:{re.escape(suffix0)}|{re.escape(suffix1)})$")
if pattern.search(self.name):
self.name = pattern.sub(rac_suffix, self.name)
self.therm = self.therm.racemize()
return self

# Private attributes
_scalers: ClassVar[Scalers] = {"therm": (lambda c, x: c * x)}

Expand Down Expand Up @@ -259,16 +271,17 @@ def fit(

# Display
def display( # noqa: PLR0913
spc: Species,
spc: Species | Sequence[Species],
*,
props: Sequence[Literal["Cv", "Cp", "S", "H", "dH"]] = ("Cp", "S", "H"),
others: Sequence[Species] = (),
others_labels: Sequence[str] = (),
T_range: tuple[float, float] = (200, 3000), # noqa: N803
units: UnitsData | None = None,
label: str = "This work",
x_label: str = "𝑇", # noqa: RUF001
label: str | Sequence[str] | None = None,
color: str | Sequence[str] | None = None,
x_label: str = "temperature",
y_labels: Sequence[str | None] | None = None,
x_unit: str | None = "K",
y_unit: str | Sequence[str | None] | None = None,
horizontal: bool = False,
) -> alt.Chart:
"""Display as an Arrhenius plot, optionally comparing to other rates.
Expand All @@ -284,15 +297,19 @@ def display( # noqa: PLR0913
:param y_label: Y-axis label
:return: Chart
"""
return spc.therm.display(
spcs = [spc] if isinstance(spc, Species) else spc
therms = [s.therm for s in spcs]
return data.display(
therm_=therms,
props=props,
others=[o.therm for o in others],
others_labels=others_labels,
T_range=T_range,
units=units,
label=label,
color=color,
x_label=x_label,
y_labels=y_labels,
x_unit=x_unit,
y_unit=y_unit,
horizontal=horizontal,
)

Expand Down
Loading
Loading