Import Simulator from library:

frompython_kinetic_simulator.simulatorimportSimulator

Set up a simulation and run it:

# initialize simulator and its temperature (in °C)sim=Simulator(T_C=50)
# add species, concentrations (M) and relative energies (kcal/mol).# Default values for energy and concentration are 0.sim.add_species("A", conc=1)
sim.add_species("B", energy=-2)
sim.add_species("B(adsorbed)", energy=-2.5)
sim.add_species("C", energy=-30)
# add reactions by specifying their activation# energy (in kcal/mol)sim.add_reaction(["A"], ["B", "B"], ts_energy=23.0)
# or their rate constant (in M^n * s^-1)sim.add_reaction(["B"], ["C"], rate=5e-4, inv_rate=0.0)
# It is also possible to enforce a K_eq directly# instead of a TS energy or rate/inverse rate.# The rate of the reaction will be assumed as fast,# to ensure equilibration at any point in time.sim.add_reaction(["B"], ["B(adsorbed)"], enforced_K_eq=2.18)
# specify the timeframe of the simulation and run itsim.run(time=10, t_units="h")
# running simulations with time=0 is also possible# to get equilibrium concentrations# sim.run(time=0)

+---+-------------------+-----------------+----------+---------+----------------------+
| # | Reaction | Faster k (s^-1) | K_eq | Pre-eq. | ΔG‡ (step, kcal/mol) |
+---+-------------------+-----------------+----------+---------+----------------------+
| 1 | B <=> B(adsorbed) | (1.00e+10) | 2.18e+00 | ✓ | (4.18) |
| 2 | A --> B + B | 1.83e-03 | 5.09e+02 | | 23.00 |
| 3 | B --> C | 5.00e-04 | 8.89e+18 | | 23.83 |
+---+-------------------+-----------------+----------+---------+----------------------+
--> Running simulation for 10 h with the LSODA method (36.0 s increments, 1000 iterations)
--> Pre-equilibrated 1 fast reactions before starting the main loop.
Iterations |##################################################| 100.0%
--> Simulation complete (0.4 s)

Show results:

sim.show()

Final Concentrations:
A : 0.00 M (0.0 % of initial conc., 100.00 % consumed)
B : 0.00 M (0.12 % total molar fraction)
B(adsorbed) : 0.01 M (0.27 % total molar fraction)
C : 1.99 M (99.61 % total molar fraction)

Only visualize the evolution of some species:

sim.show(species=("B", "B(adsorbed)"))

Final Concentrations:
B : 0.00 M (0.12 % total molar fraction)
B(adsorbed) : 0.01 M (0.27 % total molar fraction)

Simulator data is available for further manipulation:

importmatplotlib.pyplotaspltimportnumpyasnp# plot [C] vs. time in a loglog plotC_concs=sim.results["C"]
plt.loglog(sim.time_data, C_concs, label="Conc. of C over time")
plt.xlabel("ln(time (s))")
plt.ylabel("ln(conc. (M))")
# define and plot an interpolation time rangestart, end=5e1, 3e2# 1E3plt.axvspan(start, end, color="red", alpha=0.25, label="interpolation area")
# get boundary indices of arraysstart_idx=next(index+1forindex, valueinenumerate(sim.time_data[1:]) ifvalue>start)
end_idx= (
len(sim.time_data)
-2-next(indexforindex, valueinenumerate(reversed(sim.time_data[:-1])) ifvalue<end)
)
# fit line through logarithms, in the defined rangelog_time_interp=np.log(sim.time_data[start_idx:end_idx])
log_conc_interp=np.log(C_concs[start_idx:end_idx])
X=np.vstack([log_time_interp, np.ones(log_time_interp.shape[0])]).Tm, c=np.linalg.lstsq(X, log_conc_interp, rcond=None)[0]
# draw an interpolation linex_fit= (start, end)
y_fit= (start**m*np.exp(c), end**m*np.exp(c))
plt.plot(
x_fit,
y_fit,
label=f"linear fit (slope = {m:.3f})",
color="black",
linestyle="dashed",
markersize=4,
)
_=plt.legend()

It is also possible to track how much material was obtained through different pathways:

sim=Simulator()
sim.add_species("A", conc=1)
sim.add_species("B", energy=0.8)
sim.add_species("C", energy=-30)
sim.add_reaction(["A"], ["B"], ts_energy=20.0)
sim.add_reaction(["A", "B"], ["C", "C"], ts_energy=22.2, throughput_tgt="C")
sim.add_reaction(["A", "A"], ["C", "C"], ts_energy=21.5, throughput_tgt="C")
sim.run(time=24, t_units="h")
sim.show()

+---+-----------------+-----------------+----------+---------+----------------------+
| # | Reaction | Faster k (s^-1) | K_eq | Pre-eq. | ΔG‡ (step, kcal/mol) |
+---+-----------------+-----------------+----------+---------+----------------------+
| 1 | A <=> B | 5.13e-02 | 2.59e-01 | | 20.00 |
| 2 | A + B --> C + C | 1.25e-03 | 3.95e+44 | | 21.40 |
| 3 | A + A --> C + C | 1.05e-03 | 1.02e+44 | | 21.50 |
+---+-----------------+-----------------+----------+---------+----------------------+
--> Running simulation for 24 h with the LSODA method (86.4 s increments, 1000 iterations)
Iterations |##################################################| 100.0%
--> Simulation complete (0.2 s)
Final Concentrations:
A : 0.01 M (0.5 % of initial conc., 99.47 % consumed)
B : 0.00 M (0.14 % total molar fraction)
C : 0.99 M (99.34 % total molar fraction)
Reaction "A + B --> C + C" throughput is 0.228 M, 22.93 % of final "C" conc.
Reaction "A + A --> C + C" throughput is 0.766 M, 77.07 % of final "C" conc.