Reaction-rate estimation in many-body systems can be approached either by sampling stochastic trajectories or by solving the corresponding Fokker–Planck equation for the evolving probability distribution. While stochastic methods avoid explicitly representing the exponentially large configuration space and are therefore the standard classical approach, direct solution of the Fokker–Planck equation becomes intractable in high dimensions. We develop a quantum algorithm that revisits this deterministic formulation and directly estimates time-dependent reactive fluxes between prescribed reactant and product regions. After a similarity transformation, the Fokker–Planck generator becomes a negative-semidefinite self-adjoint operator admitting a sum-of-squares representation. We exploit this structure to introduce a Gaussian linear combination of Hamiltonian simulations (Gaussian-LCHS) for the non-unitary propagator, with query complexity scaling essentially as the square root of the evolution time. Directly preparing the dissipatively evolved quantum state can nevertheless incur an exponentially small postselection probability. We instead introduce a non-unitary overlap-estimation circuit that estimates the reactive flux directly as a propagator matrix element, avoiding this bottleneck. For pairwise-interacting particles, we derive explicit resource estimates and show that, under comparable worst-case analytical guarantees, the algorithm achieves an exponential improvement with particle number, a quartic improvement in accuracy, and a quadratic improvement in the time horizon relative to classical overdamped-Langevin simulation.