12. Transition State Optimization
To cite Frank Jensen: "Locating transition states (TSs) is black magic, especially in internal coordinates" and I can say this is true. The most promising TS optimizers in pysisyphus employ second derivative information (hessian) but locating TS is also possible using only first derivatives by means of the dimer method (DM).
TS optimizations should preferably be carried out in internal coordinates (type: redund|dlc). Before starting the actual TS search, one should always check the internal coordinates that pysisyphus sets up, by executing pysistrj [xyz] --internals, with [xyz] corresponding to the geometry you plan to use for the optimization. Check carefully if all supposed reaction coordinates are present. If they are missing define them manually with the add_prims key in the YAML input.
Calculation of the exact hessian can be avoided by using the DM by using rx_coords: [list of primitives] in combination with hessian_init: fischer|lindh|swart|simple. By using both options, a diagonal model hessian is calcualted and modified for use in a TS optimization.
The DM can be used with a powerful preconditioned LBFGS optimizer (type: lbfgs).
12.1. Hessian Based TS Optimization
12.1.1. YAML example
Below you can find an example YAML-input including the most important options that the user may want to modify for the Hessian-based optimizers (RS-I-RFO, RS-P-RFO, TRIM).
tsopt:
type: rsirfo|rsprfo|trim # Optimization algorithm
do_hess: True # Calculate the hessian at the final geometry
# after the optimization.
hessian_recalc: 1 # Recalculate the exact hessian every n-th cylce
# Very costly but a small number like 5 may be a good
# idea.
prim_coord: [BOND, 3, 18] # Select the mode to follow uphill by overlap with
# this primitive internal coordinate. Expects exactly
# one typed primitive.
rx_modes: [[[[BOND, 3, 18], 1]]] # Select initial mode based on overlaps with a
# constructed mode(s). Can be seen as generalization of
# prim_coord. Expects a a list of (lists of pairs, with each
# pair comprising a typed primitive and a phase factor).
#
# See examples/tsopt/{05..,06..} for examples.
#rx_coords: [[BOND, 3, 18],] # Obtain initial mode by modifying a model Hessian. Has
# to be used with 'hess_init: swart|fischer|lindh|xtb'
# to avoid calculation of exact hessian.
# Also requires 'redund' or 'dlc' to work.
#root: 0 # Follow the n-th imaginary mode uphill, minimize
# along the others.
#hessian_init: calc # Initial hessian is calculated exactly.
#hessian_update: bofill # Bofill hessian-update. Should not be modified for
# TS optimizations.
#hessian_ref: [path] # Expects a path to a precalculated hessian at any
# level of theory. The emode to maximize along is selected
# by highest overlap with imaginary mode from the reference
# hessian.
#max_micro_cycles: 50 # No. of micro cycles for the RS-variants. Does not apply
# to TRIM.
#trust_radius: 0.3 # Initial trust radius.
#trust_max: 1.0 # Max. trust radius
#trust_min: 0.1 # Min. trust radius
calc:
type: xtb
pal: 4
charge: 0
mult: 1
geom:
type: redund
fn: shaked.geom_000.xyz
add_prims: [[24, 20], ] # If using internal coordinates ALWAYS check the coordinates
# that pysisyphus generates (pysistrj [xyz] --internals). If
# some important (reaction) coordinates appears to be missing
# define them manually. In this example
# we add a bond. If additional internal coordinates can
# derived from the added primitives pysisyphus will do
# it.
Second-derivative-free TS optimization can be done using the Dimer method. Attached you can find a sample input.
Further examples for TS optimizations from .yaml input can be found here.
12.2. Dimer Method
pysisyphus implements the dimer method (DM) as calculator, wrapping another calculator
for the actual energy and gradient calculations. All DM configurations have to be done
in the calc: section. The best results are obtained by optimizing the DM with
LBFGS in connection with a preconditioner (PLBFGS), estimated from the Lindh model Hessian.
in constrast to the TS optimizers mentioned above, PBLFGS is configured in the opt
section of the YAML input. DM related information is logged to dimer.log. The current
DM orientation is saved in files with extension .N. Animated .trj files of the DM
orientation are saved to .N.trj.
12.2.1. YAML example
Below you can find an example for the DM, applied to the isomerization of HCN. Default values are commented out. See examples/tsopt/02_hcn_tsopt_dimer for the full example.
opt:
type: plbfgs # Preconditioned LBFGS
thresh: baker # Baker threshold, comparable to the Gaussian defaults
do_hess: True # Calculate Hessian at the end. May not be applicable to systems
# where Hessian calculation is infeasible.
calc:
type: dimer # Dimer calculator wrapping PySCF at HF/3-21G level of theory
calc:
type: pyscf
basis: 321g
pal: 2
charge: 0
mult: 1
#N_raw: [file] # Path to a file containing an initial dimer orientation
#length: 0.0189 # Separation of dimer images from common midpoint
#rotation_max_cycles: 15 # Maximum number of rotations
#rotation_method: fourier
#rotation_thresh: 1e-4 # rms(rotationa force) threshold
#rotation_tol: 1 # Rotation tolerance in degree. If a proposed trial
# rotation angle is below the tolerance the rotation
# will be skipped.
#rotation_max_element: 0.001 # Max. step element for "rotation_method: direct"
#rotation_interpolate: True # Interpolate force on image after rotation
#rotation_remove_trans: True # Remove overall translation from N-vector
#seed: 20182503 # Seed for the RNG for reproducability reasons.
geom:
type: cart
fn: 01_hcn.xyz
12.3. General advice for TS optimizations
Use as many exact hessians as your computational budget allows (hessian_recalc: [n])
Use tight convergence settings for your SCF. In contrast to codes like ORCA pysisyphus does not enforces this automatically
Grid-based methods (DFT, some RI-types) may need finer grids to reduce numerical noise
Try to preoptimize the TS using a cheap method (xtb) that allows hessian_recalc: 1
Decrease current & allowed trust radius (trust_radius: 0.1 and trust_max: 0.3|0.5)
Check for orrect coordinate setup
12.4. TSHessianOptimizer base class
Base class for TS optimizers that empoly a hessian.