Transition-State Search (Berny Optimizer)

The Berny optimizer supports transition-state (TS) searches via the TS task token or by setting TransitionState true in the OPT section. The implementation uses Partitioned Rational Function Optimization (P-RFO) to locate first-order saddle points.

See also: the GeomeTRIC Optimizer provides an alternative restricted-step partitioned RFO (RS-P-RFO) TS engine in translation-rotation internal coordinates. Select it with Algorithm geometric plus TransitionState true.

Theory

A transition state is a first-order saddle point on the potential energy surface (PES) — a stationary point with exactly one negative eigenvalue of the Hessian matrix. The negative eigenvalue corresponds to the reaction coordinate (the direction connecting reactant and product).

The P-RFO method partitions the Hessian eigenvectors into two subspaces:

• TS mode (mode 1, lowest eigenvalue): maximize energy along this direction
• Minimisation modes (all other modes): minimize energy along these directions

The step is computed as:

$$\Delta \mathbf{x}=\sum _i{c}_i{\mathbf{v}}_i$$

where the coefficients $c_i$ are:

$$c_i=-\frac{g_i}{\lambda -{\mu }_i}\text{for }i<\text{ModMin}[\text{TS mode: upper RFO root}]$$ $$c_i=-\frac{g_i}{\lambda -{\mu }_i}\text{for }i\ge \text{ModMin}[\text{minimisation: lower RFO root}]$$

The RFO roots are:

$${\mu }_{\text{TS}}=\frac{1}{2}\left({\lambda }_{\text{max}}+\sqrt{{\lambda }_{\text{max}}^2+4g_{\text{max}}^2}\right)[\text{upper root for TS mode}]$$ $${\mu }_{\text{min}}\text{ found by bisection to satisfy: }\sum \frac{g_i^2}{\lambda -{\mu }_i}=\mu [\text{lower root for min modes}]$$

Implementation Details

The TS search in the Berny optimizer includes:

1. Initial TS Hessian (init_ts_hessian):

   • When no analytical Hessian is provided, the Lindh model Hessian is used
   • The softest IC mode (lowest eigenvalue) is set to −0.2 a.u.
   • This ensures P-RFO has a well-defined ascent direction from the first step

2. P-RFO Step (prfo_step in step_engine.rs):

   • Eigendecomposes the Hessian: $\mathbf{H}=\mathbf{V}\mathbf{\Lambda }{\mathbf{V}}^{†}$
   • Sorts eigenvalues/eigenvectors by ascending eigenvalue
   • Projects gradient onto eigenvectors: ${\tilde{g}}_i={\mathbf{v}}_i\cdot \mathbf{g}$
   • Computes RFO roots for TS mode (upper) and minimisation modes (lower)
   • Forms step in eigenbasis, transforms back to IC space
   • Scales to trust radius if step is too large

3. Hessian Eigenvalue Check (ensure_correct_negative_eigenvalues):

   • Runs every 3 cycles during TS optimization
   • Checks that the Hessian has exactly 1 negative eigenvalue
   • If too few: flips smallest positive eigenvalues to −0.2 a.u.
   • If too many: flips least-negative eigenvalues to +0.2 a.u.

4. Bofill Hessian Update (update_hessian_bofill):

   • Blends PSB (Powell-Symmetric-Broyden) and MS (Murtagh-Sargent) updates
   • Weight: $\varphi =1-({\mathbf{s}}^{†}\cdot \mathbf{y}{)}^2/({\mathbf{s}}^{†}\cdot \mathbf{s}\cdot {\mathbf{y}}^{†}\cdot \mathbf{y})$
   • When $\varphi \to 1$: pure PSB (preserves curvature info)
   • When $\varphi \to 0$: pure MS (rank-1 correction)
   • Automatically adapts based on local surface curvature

5. TS-Specific Trust Radius:

   • Initial trust: 0.01 Å (vs 0.3 Å for minima)
   • Maximum trust: 0.03 Å (vs 0.3 Å for minima)
   • Conservative values prevent large steps near saddle points

6. TS-Specific DIIS Weighting:

   • Reaction coordinate gradient component scaled by factor 2.0
   • Minimum subspace size: 3 points (vs 2 for minima)
   • Cosine threshold: 0.5 (vs 0.0 for minima)

7. TS-Specific Output:

   • Prints Hessian negative eigenvalue count at each cycle
   • Header shows "Transition-State Optimization (Berny IC-space, P-RFO)"

TS Search Algorithm Flow

berny_optimize_ic_with_controls()
│
├─ Set TS-specific trust defaults (0.01, 0.03)
│
├─ Initialize Hessian
│   ├─ If cart_hessian provided: use it
│   └─ Else: guess_hessian_lindh()
│       └─ If transition_state: init_ts_hessian(h, 1, -0.2)
│
└─ Main loop
   │
   ├─ Compute IC gradient
   │
   ├─ TS diagnostics (print negative eigenvalue count)
   │
   ├─ Select step method
   │   ├─ If transition_state: prfo_step()
   │   └─ Else: trm_step() or gediis()
   │
   ├─ IC → Cartesian back-transform
   │
   ├─ Evaluate energy and gradient
   │
   ├─ Convergence check
   │
   ├─ Trust-radius update
   │
   ├─ Hessian update
   │   ├─ If transition_state: update_hessian_bofill()
   │   └─ Else: update_hessian_msp()
   │
   ├─ Periodic eigenvalue check (every 3 cycles)
   │   └─ ensure_correct_negative_eigenvalues(h, 1, 0.2)
   │
   └─ Accept step

Verification

A genuine TS shows exactly one imaginary frequency (printed as a negative value):

Frequencies (cm⁻¹):  -1247.3   1652.1   3825.4

Verify with a subsequent FREQ job on the converged TS geometry.