Skip to content

GeneralizedDimer: use the C++ LBFGS for the rotation step (fixes flaky test_generalized_dimer) - #48

Merged
spraharsh merged 1 commit into
masterfrom
fix-generalized-dimer
Sep 28, 2026
Merged

spraharsh merged 1 commit into
masterfrom
fix-generalized-dimer

Conversation

@spraharsh

Copy link
Copy Markdown
Collaborator

test_generalized_dimer::test1 starts from an unseeded random configuration and failed intermittently (e.g. in CI on the meson branch) with:

TypeError: LineSearchError.__init__() missing 1 required positional argument: 'value'

Cause

The dimer's rotation step (FindLowestEigenVector) minimizes the curvature along a direction with MYLBFGS(rel_energy=True), which accepts a step only if (E_new - E_old) / |E_old| <= maxErise (1e-4). When the lowest curvature approaches zero (seen at ~1e-7), ordinary finite-difference noise (~0.1% relative) exceeds that, every step is rejected, and after 10 failed line searches _lbfgs_py.py does raise LineSearchError without its required argument, hence the TypeError.

Change

  • LBFGS_CPP: get_state() / set_state() for the LBFGS memory (s, y, rho, H0, k), in C++ and Cython, with a size check
  • FindLowestEigenVector(cpp_lbfgs=True): use LBFGS_CPP and keep its memory across update_coords (as MYLBFGS did). GeneralizedDimer turns it on. findTransitionState is unchanged (still MYLBFGS); it hasn't been evaluated the same way.
  • LineSearchError is raised with a message.
  • TODO in LowestEigPot on its Python overhead: with the C++ LBFGS each curvature evaluation costs ~34 us, of which 2.5 us is the LJ18 gradient; the rest is Python/numpy/orthogopt. A C++ LowestEigPotential exists in source/pele/lowest_eig_potential.hpp but lacks rotational orthogonalization and a Cython wrapper.
  • TODO in the test with its remaining failure rate and seeds.

Evaluation

To avoid tuning on the seeds used to choose the optimizer (0-299), everything below uses held-out seeds, paired (same seeds for both rotators). A result counts as a transition state only if the exact Hessian has exactly one negative eigenvalue (rigid-body modes excluded) and rms gradient < 1e-3; "test fails" means res.success is False.

system start rotator runs verified TS (95% CI) crashes test fails
LJ18 random box (the test) MYLBFGS 2000 96.2% (95.3–97.0) 66 71
LJ18 random box (the test) LBFGS_CPP 2000 97.0% (96.2–97.7) 0 3
LJ18 minimum + noise MYLBFGS 2000 95.0% (94.0–95.9) 98 102
LJ18 minimum + noise LBFGS_CPP 2000 97.3% (96.5–97.9) 0 2

Not tuned on: LJ13, LJ38, Morse13 (rho=6) and binary LJ18, 500 seeds per start each. No crashes anywhere with LBFGS_CPP (301 → 0 over all 8000 runs), test-failure rate 334 → 67. Verified-TS rates are equal or higher within the confidence intervals, except binary LJ18 from a random box (89.8% vs 91.4%, overlapping CIs).

Known limitations, unchanged by this PR:

  • The 3 remaining LJ18 failures (seeds 10089, 10735, 11584) all lose an atom and converge to a minimum plus a free atom: the dimer has no guard against dissociation.
  • res.success only checks that the lowest curvature found is negative: 56 of the 2000 LJ18 runs "succeed" at points that are not first-order saddles (mostly two negative modes). With MYLBFGS it was 22, largely because runs that crashed now continue.

Testing

  • pytest pele/: 591 passed (built with the current setup.py), including 2 new tests: restarting LBFGS_CPP from get_state() takes the same step as the original, and set_state with the wrong size raises.

🤖 Generated with Claude Code

The dimer's rotation step minimized the curvature with MYLBFGS, which accepts a step only if the relative change in curvature is below maxErise. Near zero curvature that rejects every step, and after 10 failed line searches it raised LineSearchError without its required argument, i.e. a TypeError. test_generalized_dimer::test1 failed in 71 of 2000 held-out random starts (seeds 10000-11999); with this change it fails in 3 (an atom breaks off, see the TODO in the test).

- LBFGS_CPP: get_state/set_state for the LBFGS memory (C++ and Python)
- FindLowestEigenVector(cpp_lbfgs=True) uses LBFGS_CPP and keeps its memory across update_coords; GeneralizedDimer turns it on. findTransitionState keeps MYLBFGS.
- LineSearchError is raised with a message
- TODO on LowestEigPot's Python overhead (~34 us per call vs 2.5 us for the LJ18 gradient)

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
@codecov

codecov Bot commented Sep 27, 2026 •

Copy link
Copy Markdown

Codecov Report

❌ Patch coverage is 91.66667% with 1 line in your changes missing coverage. Please review.

Files with missing lines Patch % Lines
pele/optimize/_lbfgs_py.py 0.00% 1 Missing ⚠️

📢 Thoughts on this report? Let us know!

@spraharsh
spraharsh merged commit a6d4cdf into master Sep 28, 2026
7 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant