The model of chemical reaction networks is among the oldest and most widely studied and used in natural science. The model describes reactions among abstract chemical species, for instance $A + B \to C$, which indicates that if a molecule of type $A$ interacts with a molecule of type $B$ (the reactants), they may stick together to form a molecule of type $C$ (the product). The standard algorithm for simulating (discrete, stochastic) chemical reaction networks is the Gillespie algorithm [JPC 1977], which stochastically simulates one reaction at a time, so to simulate $\ell$ consecutive reactions, it requires total running time $Ω(\ell)$. We give the first chemical reaction network stochastic simulation algorithm that can simulate $\ell$ reactions, provably preserving the exact stochastic dynamics (sampling from precisely the same distribution as the Gillespie algorithm), yet using time provably sublinear in $\ell$. Under reasonable assumptions, our algorithm can simulate $\ell$ reactions among $n$ total molecules in time $O(\ell/\sqrt n)$ when $\ell \ge n^{5/4}$, and in time $O(\ell/n^{2/5})$ when $n \le \ell \le n^{5/4}$. Our work adapts an algorithm of Berenbrink, Hammer, Kaaser, Meyer, Penschuck, and Tran [ESA 2020] for simulating the distributed computing model known as population protocols, extending it (in a very nontrivial way) to the more general chemical reaction network setting. We provide an implementation of our algorithm as a Python package, with the core logic implemented in Rust, with remarkably fast performance in practice.

翻译：化学反应网络模型是自然科学中最古老且被广泛研究和使用的模型之一。该模型描述了抽象化学物种之间的反应，例如 $A + B \to C$，表示若一个 $A$ 型分子与一个 $B$ 型分子（反应物）相互作用，它们可能结合在一起形成一个 $C$ 型分子（产物）。模拟（离散、随机）化学反应网络的标准算法是 Gillespie 算法 [JPC 1977]，该算法每次随机模拟一个反应，因此要模拟 $\ell$ 个连续反应，其总运行时间为 $Ω(\ell)$。我们提出了首个随机模拟化学反应网络的算法，该算法能够模拟 $\ell$ 个反应，在理论上保持精确的随机动力学（从与 Gillespie 算法完全相同的分布中采样），同时使用的时间在理论上关于 $\ell$ 是亚线性的。在合理的假设下，当 $\ell \ge n^{5/4}$ 时，我们的算法可以在 $O(\ell/\sqrt n)$ 时间内模拟 $n$ 个总分子中的 $\ell$ 个反应；当 $n \le \ell \le n^{5/4}$ 时，则可在 $O(\ell/n^{2/5})$ 时间内完成。我们的工作改编了 Berenbrink、Hammer、Kaaser、Meyer、Penschuck 和 Tran [ESA 2020] 用于模拟被称为群体协议的分布式计算模型的算法，并将其（以一种非常非平凡的方式）扩展到更一般的化学反应网络场景。我们将算法实现为一个 Python 包，其核心逻辑用 Rust 编写，在实践中表现出卓越的快速性能。