About the Algorithm
The REACTER method is enabled by a novel topology-matching algorithm, referred to as the superimpose algorithm (Gissinger et al., 2017). This algorithm works behind the scenes to identify reaction sites and convert them to their post-reaction topology. Borrowing terminology from graph theory, molecules (and therefore reaction sites) can be thought of as undirected cyclic 3D graphs, and it follows that testing the equivalence of two molecules is equivalent to a graph isomorphism problem. In the case of general, arbitrary-sized graphs, the graph isomorphism problem is unsolved and even its difficulty remains unclassified in pure mathematics. However, chemical reactions generally involve no more than 30 atoms. This implies that the problem can be treated by a relatively simple, practical graph matching algorithm. These often rely on brute-force searches, but chemical intuition and other considerations result in a considerably streamlined algorithm specific to chemical structure matching. Details of the original algorithm are provided in the text and supporting information of Gissinger, Jensen and Wise (2017).
Another critical consideration of the REACTER method is how to treat the motion of 'reacting' atoms in a completely classical setting. Allowing for a reaction to occur during a running MD simulation without disrupting it is a powerful capability. It is inherently implied by such a proposal that transient unphysical behavior of atoms is unavoidable during the topology conversion. Therefore the objective of a reaction site relaxation stage is to convert reactants to products as quickly as possible while retaining simulation stability. This is additionally justified by noting that typical MD simulation time scales are much greater than those of reaction dynamics, which occur on the order of femto- or attoseconds. LAMMPS' nve/limit thermostat provides an effective option for extremely stable relaxation of reaction sites. The reaction steps tested in initial case studies require subjection to this thermostat for 30-60 time steps before being released back to the simulation's global thermostat. Additionally, fix bond/react allows for the user to easily try out custom methods of relaxing a reaction site.