Function optimize_reaction_path

Source

pub fn optimize_reaction_path(
    initial_path: &[Geometry],
    constraints: &[Constraint],
) -> Vec<Geometry>

Expand description

Optimizes a reaction path using the Nudged Elastic Band (NEB) method.

The NEB method is a powerful technique for finding minimum energy paths (MEPs) between known reactant and product structures. It optimizes a chain of intermediate images (geometries) connected by springs to find the most probable reaction pathway.

§Theoretical Foundation

§The NEB Method

NEB balances two competing forces on each intermediate image:

Spring Forces (F^spring): Keep images evenly distributed along the path
True Forces (F^true): Drive images toward lower energy regions

The total NEB force on image i is:

F_i^NEB = F_i^spring + F_i^⊥

where F_i^⊥ is the component of the true force perpendicular to the path.

§Spring Force Calculation

The spring force along the path tangent τ̂ is:

F_i^spring = k[(|R_{i+1} - R_i| - |R_i - R_{i-1}|)] τ̂_i

where:

k is the spring constant (typically 0.1-1.0 eV/Angstrom²)
R_i represents the coordinates of image i
τ̂_i is the normalized tangent vector at image i

§True Force Projection

The perpendicular component of the true force is:

F_i^⊥ = F_i^true - (F_i^true · τ̂_i) τ̂_i

This ensures forces don’t interfere with spring-controlled spacing.

§Algorithm Overview

Initialization: Start with initial path (linear interpolation or guess)
Force Calculation: For each intermediate image:
- Calculate tangent vector from neighboring images
- Compute spring forces along the tangent
- Project true forces perpendicular to path
Geometry Update: Move each image according to total NEB force
Convergence Check: Continue until forces are below threshold
Constraint Application: Apply any geometric constraints if specified

§Arguments

initial_path - Starting path as a series of geometries (reactant → product)
constraints - Optional geometric constraints to maintain during optimization

§Returns

Returns an optimized Vec<Geometry> representing the minimum energy path. The first and last geometries (endpoints) remain fixed during optimization.

§Implementation Details

§Current Implementation (Simplified)

This implementation uses a simplified NEB approach:

Spring forces only: True energy gradients not included
Fixed endpoints: Reactant and product geometries unchanged
Steepest descent: Simple optimization without advanced methods
Basic constraints: Limited constraint handling

§Parameters

Spring constant: k = 0.1 (relatively soft springs)
Step size: α = 0.01 (conservative for stability)
Max iterations: 100 (sufficient for most systems)
Convergence threshold: 1e-3 (moderate precision)

§Applications

Transition State Location: Find saddle points along reaction paths
Reaction Mechanism Studies: Understand detailed reaction pathways
Activation Energy Calculation: Determine energy barriers
Conformational Transitions: Study large-scale molecular motions

§Examples

use omecp::reaction_path::optimize_reaction_path;

// Create initial path by linear interpolation
let initial_path = vec![reactant_geom, intermediate_geom, product_geom];

// Optimize the path using NEB
let optimized_path = optimize_reaction_path(&initial_path, &[]);

println!("Optimized path has {} images", optimized_path.len());

§Advanced NEB Variants (Not Implemented)

Climbing Image NEB (CI-NEB): Drives highest energy image to saddle point
Adaptive NEB: Automatically adjusts number of images
String Method: Alternative path optimization approach
Growing String Method: Dynamically extends path length

§References

Henkelman, G.; Jónsson, H. J. Chem. Phys. 2000, 113, 9978-9985.
Henkelman, G.; Uberuaga, B. P.; Jónsson, H. J. Chem. Phys. 2000, 113, 9901-9904.

optimize_reaction_path

Function optimize_reaction_path Copy item path