 codex/esipt-population-rise.py

"""Illustrative irreversible excited-state proton-transfer kinetics (Python 3.14)."""
from pathlib import Path
import numpy as np
import matplotlib
matplotlib.use("Agg")
import matplotlib.pyplot as plt


def populations(t, a=10.0, b=1.0, transfer=8.0):
    normal = np.exp(-a * t)
    if np.isclose(a, b):
        tautomer = transfer * t * normal
    else:
        tautomer = transfer * (normal - np.exp(-b * t)) / (b - a)
    return normal, tautomer


def main():
    t = np.linspace(0.0, 3.0, 1000)
    normal, tautomer = populations(t)
    assert np.all(normal >= 0) and np.all(tautomer >= 0)
    assert np.all(normal + tautomer <= 1.0 + 1e-12)
    peak_time = np.log(10.0) / 9.0
    peak = populations(np.array([peak_time]))[1][0]
    fig, ax = plt.subplots(figsize=(8.0, 4.2), dpi=100, facecolor="white")
    ax.set_facecolor("white")
    ax.plot(t, normal, color="#235d9c", lw=2.3, label=r"Normal state $N^*$")
    ax.plot(t, tautomer, color="#a64919", lw=2.3, label=r"Transferred state $T^*$")
    ax.plot(peak_time, peak, "o", color="#a64919", ms=5)
    ax.annotate(f"Tautomer peaks at {peak_time:.3f} ns", xy=(peak_time, peak),
                xytext=(0.85, 0.71), arrowprops={"arrowstyle": "->", "color": "#555555"},
                fontsize=10)
    ax.set(xlim=(0, 3), ylim=(0, 1.04), xlabel="Time after ideal light pulse (ns)",
           ylabel="Population / initially excited population")
    ax.set_title("An excited-state tautomer can rise after a light pulse", fontsize=13, pad=26)
    ax.text(0.5, 1.025, r"Illustrative model: $a=10$, $b=1$, $k_p=8$ ns$^{-1}$; no reverse transfer",
            transform=ax.transAxes, ha="center", va="bottom", fontsize=10, color="#444444")
    ax.grid(alpha=0.20)
    ax.legend(loc="upper right", frameon=False)
    fig.tight_layout()
    fig.savefig(Path(__file__).with_suffix(".png").name, facecolor="white", dpi=100)
    plt.close(fig)


if __name__ == "__main__":
    main()