Management of groundwater resources and remediation of groundwater pollution require reliable quantification of contaminant dy-namics in natural aquifers, which can involve complex chemical dynamics and challenge traditional modeling approaches. The kinetics of chemical reactions in groundwater are well known to be controlled by medium heterogeneity and reactant mixing, motivating the development of particle-based Lagrangian approaches. Previous Lagrangian solvers have been limited to fundamental bimolecular reactions in typically one-dimensional porous media. In contrast to other existing studies, this study developed a fully Lagrangian framework, which was used to simulate diffusion-controlled, multi-step reactions in one-, two-, and three-dimensional porous media. The interaction radius of a reactant molecule, which controls the probability of reaction, was derived by the agent-based approach for both irreversible and reversible reactions. A flexible particle tracking scheme was then developed to build trajectories for particles undergoing mixing-limited, multi-step reactions. The simulated particle dynamics were checked against the kinetics for diffusion-controlled reactions and thermodynamic well-mixed reactions in one- and two-dimensional domains. Applicability of the novel simulator was further tested by (1) simulating precip-itation of calcium carbonate minerals in a two-dimensional medium, and (2) quantifying multi-step chemical reactions observed in the laboratory. The flexibility of the Lagrangian simulator allows further refinement to capture complex transport affecting chemical mixing and hence reactions.