← All findings

F0002: mepd channels finds g-xTB transition states below the reference ones in Reaction-QM

Status: verified · Version 5, updated 2026-10-24 · Jan Estrada Pabon (ORCID) · 2026-10-24 · no DOI yet

Finding. Searching with mepd channels on the g-xTB surface, 8,737 of the 31,617 Reaction-QM reactions of the F0001 benchmark set gave a verified route, every step of it reaching its minima by its own IRC, whose highest transition state lies more than 1 kcal/mol below the re-converged reference one (27.6%; median 5.04 kcal/mol lower). For 5,129 (16.2%) no route was found.

Updated in version 5 (2026-10-24): corrected. Route steps must now reach their minima by their own IRC, and the benchmark set (F0001, version 3) retries a reference IRC that stops short before leaving a reaction out. 8,737 of the 31,617 reactions now in the set have a lower route (version 4: 8,150 of 29,116); in 315 the best route changed (159 now have none). Version 4's counts should not be used. How this Finding changed

Background

mepd is a tool for mapping reaction mechanisms. Given a reactant and a product, its search, mepd channels, looks for the transition state between them: the top of the energy hill the reaction has to climb. A lower hill means a faster reaction, so the lowest transition state connecting the two is the one that matters.

We test the search on reactions whose transition state is already known, from Reaction-QM. The search runs on g-xTB, a fast approximate quantum-chemistry method, so each known transition state was first re-converged on g-xTB and kept only when it still leads to both of the dataset's endpoints, by its own IRC (Finding F0001, version 3). The question here: does the search find a route that is lower than the known one?

Approach

For each reaction of the F0001 benchmark set, mepd channels built several starting arrangements of the reactant and product and searched for paths between them on g-xTB.¹ Every transition state it reports was checked: optimized, then followed downhill to confirm it joins the requested reactant and product, either directly or through intermediates.² Each step's IRC must reach its minima itself: when minimizing an IRC end turned it into another molecule, the IRC was followed again with a ten times stricter stopping point, and that IRC decides.³

We compared the energy of the best route's highest transition state with that of the known transition state, both on g-xTB. A route counts as lower when it is more than 1 kcal/mol below.⁴ All 31,617 reactions of the set were searched.⁵ Version 4 had neither the IRC retries nor the step check.⁶

Results

The 10 single-step and 10 multistep routes furthest below the known transition state, beside the dataset's own paths:

Most reactions with no route are a limit of the search budget, not evidence that no route exists: the pair of the reference IRC's own two ends was usually not searched in time.⁹ These are counts of what one search budget found, so they are a floor, not a ceiling.¹⁰ Whether the lower routes stay lower at a higher level of theory was not tested here.

History

How this Finding has changed as more data came in, newest first.

Errors corrected in version 5:

• git commit 87d9f69 of mepd-workflows: routes were counted although a step's IRC reached its intermediate or end only after minimization turned the IRC end into another molecule
• git commit 87d9f69 of mepd-workflows: reactions whose reference IRC stopped short at fmax 0.1 eV/Å, so that minimizing an end changed its species, were left out without retrying a stricter IRC

Errors corrected in version 3:

• git commit 5f3b655 of mepd-workflows: the benchmark set kept reactions whose reference IRC reached the dataset's endpoints only after minimization changed an IRC end into another species

Errors corrected in version 2:

• git commit 5b18bda of mepd: the path search passed coordinates to geodesic interpolation in the wrong unit (bohr for ångström)
• git commit 0594999 of mepd: reference transition states and IRC ends were not checked strictly (convergence, frequencies, minimized ends, stereochemistry)

Version 4 (2026-10-22) said: Outdated: corrected in version 5. Read version 4, with its data.

Searching with mepd channels on the g-xTB surface, 8,150 of the 29,116 Reaction-QM reactions of the F0001 benchmark set gave a verified route whose highest transition state lies more than 1 kcal/mol below the re-converged reference one (28.0%; median 4.98 kcal/mol lower). For 4,459 (15.3%) no route was found.

Version 3 (2026-10-22) said: Outdated: corrected in version 5. Read version 3, with its data.

Searching with mepd channels on the g-xTB surface, 8,150 of the 29,116 Reaction-QM reactions of the F0001 benchmark set gave a verified route whose highest transition state lies more than 1 kcal/mol below the re-converged reference one (28.0%; median 4.98 kcal/mol lower). For 4,459 (15.3%) no route was found.

Version 2 (2026-10-22) said: Outdated: corrected in version 3. Read version 2, with its data.

Searching with mepd channels on the g-xTB surface, 8,851 of the 31,914 Reaction-QM reactions of the F0001 benchmark set gave a verified route whose highest transition state lies more than 1 kcal/mol below the re-converged reference one (27.7%; median 5.11 kcal/mol lower). For 5,153 (16.1%) no route was found.

Version 1 (2026-10-05) said: Outdated: corrected in version 2. Read version 1, with its data.

Searching with mepd channels on the g-xTB surface found a verified route whose highest transition state lies more than 1 kcal/mol below the dataset's own transition state (re-converged at g-xTB) for 103 of 500 Reaction-QM reactions (median 6.02 kcal/mol lower) and for 12 of 66 Geodesic-set reactions (median 15.72 kcal/mol lower).

Details

1. The search

mepd channels, with the settings of the original 500-reaction run. Each setting is explained in the mepd manual's page on channel searches (online; a copy is in this Finding's download as reproduce/mepd-manual-channels.html):

• 5 conformers per side (n_conformers): up to five RDKit conformers each of the reactant and of the product, chosen to differ in shape, each minimized on g-xTB. Every reactant conformer is paired with every product conformer.
• 5 path searches per mechanism (pairs_per_mechanism): the pairs are grouped by mechanism (which bonds break and which form, from SLAPMapper's atom mapping), and in each group only the five pairs whose two ends are closest in shape (endpoint-rmsd) get a path search.
• Direct paths only (direct_only): when a path turns out to pass through intermediates, only the pieces that start or end at the requested reactant or product are searched further.
• 30-minute search budget (search_budget): after 30 minutes no new path search starts and paths are no longer split; the transition-state optimizations and IRCs of what was found still run.
• 2-hour hard limit (hard time limit): the benchmark script stops a reaction's search after 2 hours, and the reaction then has no result.
• GSM with 20 nodes and at most 80 iterations (GSM settings) as the path method, one search per core.
• Checking intermediates (hessian_minimum_frequency_cutoff, hessian_minima_rescue_displacement): a path is split at an apparent intermediate only if it has no frequency below −20 cm⁻¹; one that fails is pushed 0.3 bohr along its lowest mode, each way, and minimized again.

The reactant and product are the two ends of the re-converged reference transition state's IRC, each minimized (F0001), so every search starts from the structures the reference is compared on. The searches ran with mepd git commits 0594999 (the first 1,248), 9aee085 and da31dcb (the rest, among them all reactions new in version 5); the later two only use less memory and give the same routes and energies.

2. What counts as a found route

Only routes mepd itself verified count. A single-step route is one transition state (exactly one imaginary frequency) whose IRC, with both ends minimized, reconnects the requested reactant and product species (bond graph and stereochemistry). A multistep route is a chain of such transition states through intermediates from reactant to product; its height is its highest transition state. Routes that lead to other species are not counted.

3. Each step must reach its minima by its own IRC

mepd channels files each step under the minima its IRC ends minimize to, without checking that minimizing kept each end the same molecule. We checked every step of every route found, lowest route first, with the rule of the benchmark set: when minimizing an IRC end changed it into another molecule, the step's IRC was run again from its transition state with fmax 0.01 eV/Å (at most 3000 steps). The step stands when both ends of that IRC keep their species on minimization and are the same two molecules as before. A route counts only when all its steps stand, and the best route is the lowest that counts. 2,558 step IRCs in 2,361 reactions were retried. The best route changed in 315 reactions; in 159 no route was left, and they count as no route found. routes_dropped and step_ircs_retried in data/search_outcomes.csv give the counts per reaction.

4. Comparing energies

Both barriers are measured from one floor per reaction, the lowest g-xTB energy of the reactant and its search conformers, so the amount lower is the reference transition state's energy minus that of the best route's highest one. Outcomes: lower (more than 1 kcal/mol below), about the same (within 1), only higher (every route more than 1 kcal/mol above). An IRC end of a route found can lie slightly below the floor (another arrangement of the reactant molecules); the barriers shown are then a little high, but the amounts lower do not depend on the floor.

5. Which reactions

The F0001 benchmark set: of Reaction-QM's B3LYP-RXN reactions, the 31,617 whose reference transition state, re-converged on g-xTB, leads by its IRC to both of the dataset's endpoints, with the same bond graph and stereochemistry. The IRC itself has to get there: each IRC end is minimized, and the reaction is kept only when that minimization leaves it the same species. Where it did not, the IRC was followed again with fmax 0.01 eV/Å instead of 0.1 (at most 3000 steps per direction), and that IRC decides (F0001, version 3). In 337 more reactions, even with the stricter IRC, the minimized ends matched the dataset's only because minimizing changed an IRC end into another molecule. Mostly the IRC stopped far from a minimum: in 213 an end dropped more than 20 kcal/mol when minimized. Some sit on very soft saddles (56 with an imaginary mode above −200 cm⁻¹; RXN_0000011948: −51 cm⁻¹, both ends dropping about 130 kcal/mol). Then the reference transition state is not the reaction's transition state, so these reactions are left out; version 2 of this Finding counted them. The rule is cautious: in 65 of them the saddle is clear and no end drops more than 5 kcal/mol, which may still be the same reaction; they are left out too. All are listed in data/left_out_irc.csv. All reactions in the set are neutral closed-shell singlets of 5 to 33 atoms.

6. What version 4 got wrong

Version 4 counted routes with a step whose IRC reached its intermediate or end only after minimization turned the IRC end into another molecule, so the minimization, not the transition state, did part of the chemistry. Its benchmark set also left out every reaction whose reference IRC had that problem, without trying a stricter IRC (F0001, version 3, corrects this). Version 4's counts (29,116 reactions in the set) should not be used. Its searches are reused unchanged where the reaction's reference is unchanged. The reactions the stricter IRC brought into the set were searched with the same settings; where such a reaction had been searched for version 2, from the ends of the first IRC, that search was set aside. Five reactions whose search files were lost after their outcome was recorded in version 4 were searched again.

7. No route found

No transition state connecting the requested reactant and product was found within the budget: either no path search led to one, every transition state found led elsewhere, or every route found had a step whose IRC did not reach its minima by itself.³ It is not evidence that no such route exists. These reactions are listed in data/search_outcomes.csv with status no_route.

8. Searches run twice

Eight searches first stopped early: three were killed for running out of memory, three reached the 2-hour limit and two stopped with an error. All eight were run again with the same settings, with mepd git commit da31dcb, which only lowers a search's memory use. Six finished. The other two stopped with the same error again: every path search between their reactant and product failed, so mepd had no result to compare and reports an error instead of "no route found". Of the searches new in version 5, 1 stopped with the same error. They are listed in data/search_outcomes.csv with status error.

9. Why no route was found

Each search builds conformers of the two ends, pairs them, and runs a path search pair by pair until its 30-minute budget is used up. One candidate pair is the reference IRC's own two ends, between which the reference transition state is known to lie. For the 4,970 reactions where no route was found at all (leaving out the 159 whose routes were all set aside), that pair:

A search keeps at most five pairs per bond-change mechanism, ranked by how well the two ends map onto each other, so the pair of IRC ends can be left out; in 261 of these reactions an IRC end was not among the conformers the search kept at all. 4,624 of the 4,970 searches (93.0%) used up the budget, and 4,810 found transition states that lead to other species. So most "no route" outcomes mean the right candidate was buried behind the budget, not that the method cannot find these transition states; how to order and budget the pairs is open work. Pairs are compared by geometry (sorted interatomic distances, 0.02 Å), since a search renumbers atoms when it maps the two ends. Per reaction: data/search_diagnostics.csv.

10. Where this holds

Level of theory. g-xTB 2.0.1 throughout: path search (GSM), TS optimization and IRC (Sella through ASE and mepd's g-xTB engine); gas phase. Reference TSs are Reaction-QM's B3LYP TSs re-converged at g-xTB (F0001).

Software.

• mepd: git commits 0594999, 9aee085 and da31dcb (same results; the later two use less memory)
• g-xTB: 2.0.1 (binary gxtb-xtb)
• GSM: patched molecularGSM (fork jandestrada/molecularGSM, branch mepd, git commit b21a100)

Chemical domain.

• Dataset: Reaction-QM B3LYP-RXN: the F0001 benchmark set (31,617 reactions)
• Elements: H, C, N, O, F, Si, P, S, Cl
• Total charges: 0
• Spin multiplicities: 1
• Atoms: 5–33
• Notes: All neutral closed-shell singlets.

Settings the result depends on.

• Search: mepd channels: 5 RDKit conformers per side, direct paths only, 5 pairs per mechanism, RMSD cutoff 0.1, 1800 s search budget, 7200 s hard limit, Hessian minimum check (−20 cm⁻¹, 0.3 bohr rescue); GSM 20 nodes, 80 iterations
• Atom mapping: SLAPMapper, endpoint-RMSD metric
• Endpoints: the minimized IRC ends of the re-converged reference TS (F0001 version 3: IRC retried at fmax 0.01 eV/Å, 3000 steps, where minimizing changed an end's species)
• Route steps: a route counts only when every step's IRC ends keep their species on minimization (a failing step's IRC retried at fmax 0.01 eV/Å must then reach the same two species)
• Lower: best route's highest TS more than 1 kcal/mol below the reference TS, both from one shared floor

Tested on 31,617. All 31,617 reactions of the F0001 benchmark set, version 3 (reference IRC, retried stricter where needed, reaches both dataset endpoints without minimization changing an end's species).

Not tested.

• Whether the lower transition states are also lower at a higher level of theory.
• Reactions where the search found no route: a lower transition state may exist there too.
• Longer searches, more conformers or other path methods; the counts are what this budget found.
• Whether the reactions with no route have one that a larger budget, or searching the pair of IRC ends first, would find.

Known exceptions.

• An IRC end of a route found can lie slightly below the floor; barriers are then a little high, amounts lower are unaffected.
• 3 reactions stop with an error: all their path searches fail, so they have no result.
• Reactions whose reference IRC reached the dataset's endpoints only after minimization changed an end's species, even with the stricter IRC (337, listed in data/left_out_irc.csv), are not in the set.
• The step check compares bond graphs and stereochemistry; a few steps fail it only because perception is undecided at a distorted IRC end (a bond at the length threshold, bond orders not assigned), and those routes are dropped.

Reproducibility: statistical. Rerunning the searches is not expected to find exactly the same routes: they depend on random conformers, timing and the code version.

11. Data

data/benchmark_set.csv columns:

• reaction_id: Reaction-QM id
• n_atoms: atoms

data/search_outcomes.csv columns:

• reaction_id: Reaction-QM id
• n_atoms: atoms
• status: done (a route counts) | no_route (none found, or every route found had a step that failed the IRC check) | timeout (stopped at the 2-hour limit) | error (every path search failed, so the search had no result to compare)
• outcome: lower (> 1 kcal/mol below the reference) | matched (within 1) | higher (> 1 above); empty without a route
• best_kind: direct (one step) | multistep (its highest TS counts)
• ref_barrier_kcal: re-converged reference TS above the floor, kcal/mol
• best_barrier_kcal: best route's highest TS above the floor, kcal/mol
• lower_by_kcal: ref_barrier_kcal minus best_barrier_kcal (positive = found lower)
• n_direct: single-step routes found
• n_multistep: multistep routes found
• seconds: search wall time, s
• dataset_barrier_kcal: Reaction-QM's own B3LYP barrier, kcal/mol (for reference only)
• screen_barrier_kcal: g-xTB barrier of the re-converged reference TS from its own IRC end (F0001), kcal/mol
• routes_dropped: routes found below the best one that counts, set aside because a step's IRC did not reach its minima by itself
• step_ircs_retried: route-step IRCs retried at fmax 0.01 eV/Å for this reaction

data/left_out_irc.csv columns:

• reaction_id: Reaction-QM id
• n_atoms: atoms
• ts_imaginary_cm: imaginary frequency of the re-converged reference TS, cm⁻¹
• start_end_drop_kcal: energy drop when the reactant-side IRC end was minimized, kcal/mol
• end_end_drop_kcal: the same for the product-side end, kcal/mol
• start_end_changed_species: minimizing the reactant-side IRC end changed its species
• end_end_changed_species: minimizing the product-side IRC end changed its species
• v2_status: the search's status in version 2
• v2_outcome: its outcome in version 2
• v2_lower_by_kcal: version 2's amount lower, kcal/mol (against a reference that is not this reaction's TS)

data/search_diagnostics.csv columns:

• reaction_id: Reaction-QM id
• status: done (a route was found) | no_route
• readable: 1 if the search folder could be read
• n_pairs_possible: conformer pairs the search could form
• n_pairs_queued: pairs it chose to search
• irc_ends_in_pools: 1 if both IRC ends were among the conformers kept
• irc_pair_queued: 1 if the pair of IRC ends was chosen
• irc_pair_path_done: 1 if its path search finished
• irc_pair_path_split: 1 if that path went through intermediates
• budget_reached: 1 if the search used up its 30-minute budget
• legs_not_run: path legs through intermediates not followed (direct-only searches)
• n_offtarget: transition states found that lead to other species

check.py recomputes all 36 numbers in this report from these files and compares them with expected_numbers.json.

12. Reproducing it

Software. Python 3.12, mepd, g-xTB 2.0.1, patched molecularGSM, RDKit, Sella; one core per search.

mepd version. Git commit da31dcb, the start of the commit hash that identifies this exact version of the mepd code (git checkout da31dcb).

Cost. About 18,000 core-hours for the searches (median 37 minutes per reaction on one core); the route-step check needs a few core-minutes per reaction where a step's IRC is retried.

Download. F0002-v5.zip holds the data, the check script, the rerun scripts and the input structures.

Check the numbers (a few seconds):

python check.py

Rerun the searches. This needs uv and access to mepd (git commit da31dcb; mepd is not public yet, so uv fetches it over HTTPS with your GitHub credentials, e.g. after gh auth setup-git). It also needs g-xTB 2.0.1 on PATH as gxtb-xtb, and the patched GSM from jandestrada/molecularGSM (branch mepd) on PATH as gsm. Each search takes up to 2 hours on one core (median 37 minutes), so all 31,617 took about 18,000 core-hours; --jobs sets how many run at once. rescore then checks every step of every route found and retries the IRC of a step whose end changed species (a few core-minutes per reaction that needs it).

cd reproduce/rerun
tar xJf ../../data/large/inputs.tar.xz        # -> ws/
export OMP_NUM_THREADS=1
uv run --no-config benchmark.py search ws --profile mepd_profile.toml --options search_options.toml \
  --timeout 7200 --jobs 8
uv run --no-config benchmark.py rescore ws --timeout 3600 --jobs 8
python outcomes.py ws rerun_search_outcomes.csv

Then compare rerun_search_outcomes.csv with data/search_outcomes.csv. A run can be stopped and started again: reactions already searched are skipped. To try a few reactions first, keep only their lines in ws/screen/results.jsonl before the search. The routes found may not match exactly: they depend on random conformers and timing. (benchmark.py can also send the searches to the authors' cluster, --backend afh; that is not needed here.)

Reproductions.

• Passed on 2026-10-24, rerun by Jan Estrada Pabon. quick rerun exit 0

13. Provenance

Captured: recorded by the run itself (provenance.json). Each search stage records its mepd provenance (code version, settings, programs). The retried route-step IRCs ran with mepd git commit da31dcb; the comparisons of the route-step check ran with mepd git commit 257dbb7, whose changes since da31dcb concern metal centres only (none in this set).

14. Sources

Data this builds on.

• F0001: g-xTB benchmark set of Reaction-QM transition states
• Reaction-QM data: doi:10.5281/zenodo.18551029, license CC-BY-4.0

Cited.

• Reaction-QM (B3LYP-RXN subset), Sci Data 2026: doi:10.1038/s41597-026-07325-w

Software used.

• g-xTB: https://github.com/grimme-lab/g-xtb

Other Findings.

• Builds on F0001