17.1.1.10. pysisyphus.irc package

17.1.1.10.1. Submodules

17.1.1.10.2. pysisyphus.irc.DWI module

class pysisyphus.irc.DWI.DWI(n=4, maxlen=2)[source]

Bases: object

__init__(n=4, maxlen=2)[source]

Distance weighted interpolation.

dump(fn)[source]
static from_h5(fn)[source]
interpolate(at_coords, gradient=False)[source]

See [1] Eq. (25) - (29)

update(coords, energy, gradient, hessian)[source]
pysisyphus.irc.DWI.taylor(energy, gradient, hessian, step)[source]

Taylor series expansion of the energy to second order.

pysisyphus.irc.DWI.taylor_grad(gradient, hessian, step)[source]

Gradient of a Taylor series expansion of the energy to second order.

17.1.1.10.3. pysisyphus.irc.DampedVelocityVerlet module

class pysisyphus.irc.DampedVelocityVerlet.DampedVelocityVerlet(geometry, v0=0.04, dt0=0.5, error_tol=0.003, max_cycles=150, **kwargs)[source]

Bases: IRC

damp_velocity(velocity)[source]
estimate_error(new_mw_coords)[source]
mw_grad_to_acc(mw_grad)[source]

Takes care of the units for the mass-weighted gradient. Converts units of a mass-weighted gradient [Hartree/(Bohr*amu)] to units of acceleration [sqrt(amu)*Bohr/fs²]. The 1e30 comes from converting second² to femto second².

prepare(direction)[source]
step()[source]

17.1.1.10.4. pysisyphus.irc.Euler module

class pysisyphus.irc.Euler.Euler(geometry, step_length=0.01, **kwargs)[source]

Bases: IRC

step()[source]

17.1.1.10.5. pysisyphus.irc.EulerPC module

class pysisyphus.irc.EulerPC.EulerPC(*args, hessian_recalc=None, hessian_update='bofill', hessian_init='calc', max_pred_steps=500, loose_cycles=3, dump_dwi=False, scipy_method=None, corr_func='mbs', **kwargs)[source]

Bases: IRC

corrector_step(init_mw_coords, step_length, dwi)[source]
get_integration_length_func(init_mw_coords)[source]
prepare(*args, **kwargs)[source]
scipy_corrector_step(init_mw_coords, step_length, dwi)[source]

Solve IRC equation dx/ds = -g/|g| on DWI PES in mass-weighted coordinates. Integration done until self.step_length in unweighted coordinates is achieved.

step()[source]

17.1.1.10.6. pysisyphus.irc.GonzalezSchlegel module

class pysisyphus.irc.GonzalezSchlegel.GonzalezSchlegel(geometry, max_micro_cycles=20, micro_step_thresh=0.001, hessian_recalc=None, line_search=False, **kwargs)[source]

Bases: IRC

micro_step(counter)[source]

Constrained optimization on a hypersphere.

perp_component(vec, perp_to)[source]
postprocess()[source]
step()[source]

17.1.1.10.7. pysisyphus.irc.IMKMod module

class pysisyphus.irc.IMKMod.IMKMod(geometry, corr_first=True, corr_first_size=0.5, corr_first_energy=True, corr_bisec_size=0.0025, corr_bisec_energy=True, **kwargs)[source]

Bases: IRC

fit_grad_and_energies(energy0, grad0, energy1, direction)[source]
fit_parabola(x, y)[source]
step()[source]

17.1.1.10.8. pysisyphus.irc.IRC module

class pysisyphus.irc.IRC.IRC(geometry, step_length=0.1, max_cycles=125, downhill=False, forward=True, backward=True, root=0, hessian_init=None, displ='energy', displ_energy=0.001, displ_length=0.1, displ_third_h5=None, rms_grad_thresh=0.001, hard_rms_grad_thresh=None, energy_thresh=1e-06, imag_below=0.0, force_inflection=True, check_bonds=True, out_dir='.', prefix='', dump_fn='irc_data.h5', dump_every=5)[source]

Bases: object

__init__(geometry, step_length=0.1, max_cycles=125, downhill=False, forward=True, backward=True, root=0, hessian_init=None, displ='energy', displ_energy=0.001, displ_length=0.1, displ_third_h5=None, rms_grad_thresh=0.001, hard_rms_grad_thresh=None, energy_thresh=1e-06, imag_below=0.0, force_inflection=True, check_bonds=True, out_dir='.', prefix='', dump_fn='irc_data.h5', dump_every=5)[source]

Base class for IRC calculations.

Parameters:
  • geometry (Geometry) -- Transtion state geometry, or initial geometry for downhill run.

  • step_length (float, optional) -- Step length in unweighted coordinates.

  • max_cycles (int, optional) -- Positive integer, controlloing the maximum number of IRC steps taken in a direction (forward/backward/downhill).

  • downhill (bool, default=False) -- Downhill run from a non-stationary point with non-vanishing gradient. Disables forward and backward runs.

  • forward (bool, default=True) -- Integrate IRC in positive s direction.

  • backward (bool, default=True) -- Integrate IRC in negative s direction.

  • root (int, default=0) -- Use n-th root for initial displacement from TS.

  • hessian_init (str, default=None) -- Path to Hessian HDF5 file, e.g., from a previous TS calculation.

  • displ (str, one of ("energy", "length", "energy_cubic")) -- Controlls initial displacement from the TS. 'energy' assumes a quadratic model, from which a step length for a given energy lowering (see 'displ_energy') is determined. 'length' corresponds to a displacement along the transition vector. 'energy_cubic' considers 3rd derivatives of the energy along the transition vector.

  • displ_energy (float, default=1e-3) -- Required energy lowering from the TS in au (Hartree). Used with 'displ: energy|energy_cubic'.

  • displ_length (float, default=0.1) -- Step length along the transition vector. Used only with 'displ: length'.

  • displ_third_h5 (str, optional) -- Path to HDF5 file containing 3rd derivative information. Used with 'displ: energy_cubic'.

  • rms_grad_thresh (float, default=1e-3,) -- Convergence is signalled when to root mean square of the unweighted gradient is less than or equal to this value

  • energy_thresh (float, default=1e-6,) -- Signal convergence when the energy difference between two points is equal to or less than 'energy_thresh'.

  • imag_below (float, default=0.0) -- Require the wavenumber of the imaginary mode to be below the given threshold. If given, it should be a negative number.

  • force_inflection (bool, optional) -- Don't indicate convergence before passing an inflection point.

  • check_bonds (bool, optional, default=True) -- Report whether bonds are formed/broken along the IRC, w.r.t the TS.

  • out_dir (str, optional) -- Dump everything into 'out_dir' directory instead of the CWD.

  • prefix (str, optional) -- Short string that is prepended to all files that are created by this class, e.g., trajectories and HDF5 dumps.

  • dump_fn (str, optional) -- Base name for the HDF5 files.

  • dump_every (int, optional) -- Dump to HDF5 every n-th cycle.

property coords
dump_data(dump_fn=None, full=False)[source]
dump_ends(path, prefix, coords=None, trj=False)[source]
property energy
get_conv_fact(mw_grad, min_fact=2.0)[source]
get_endpoint_and_ts_geoms()[source]
get_full_irc_data()[source]
get_irc_data()[source]
get_path_for_fn(fn)[source]
property gradient
initial_displacement()[source]

Returns non-mass-weighted steps in +s and -s direction for initial displacement from the TS. Earlier version only returned one step, that was later multiplied by either 1 or -1, depending on the desired IRC direction (forward/backward). The current implementation directly returns two steps for forward and backward direction. Whereas for plus and minus steps for displ 'length' and displ 'energy'

step_plus = -step_minus

is valid, it is not valid for dipsl 'energy_cubic' anymore. The latter step is formed as

x(ds) = ds * v0 + ds**2 * v1

so

x(ds) != -x(ds)

as

ds * v0 + ds**2 * v1 != -ds * v0 - ds**2 * v1 .

So, all required step are formed directly and later used as appropriate.

See

https://aip.scitation.org/doi/pdf/10.1063/1.454172 https://pubs.acs.org/doi/10.1021/j100338a027 https://aip.scitation.org/doi/pdf/10.1063/1.459634

irc(direction)[source]
log(msg)[source]
property m_sqrt
mass_weigh_hessian(hessian)[source]
property mw_coords
property mw_gradient
postprocess()[source]
prepare(direction)[source]
report_bonds(prefix, bonds)[source]
report_conv_thresholds()[source]
run()[source]
set_data(prefix)[source]
unweight_mw_grad(vec)[source]
valid_displs = ('energy', 'length', 'energy_cubic')

17.1.1.10.9. pysisyphus.irc.IRCDummy module

class pysisyphus.irc.IRCDummy.IRCDummy(all_coords, atoms, forward=True, backward=True, downhill=False)[source]

Bases: object

all_coords: list
atoms: tuple
backward: bool = True
downhill: bool = False
forward: bool = True

17.1.1.10.10. pysisyphus.irc.Instanton module

class pysisyphus.irc.Instanton.Instanton(images, calc_getter, T)[source]

Bases: object

property P
P_bh

the first image is connected to the last image. At k=0 the index k-1 will be -1, which points to the last image.

Below we pre-calculate some indices (assuming N images).

unshifted indices ks: k = {0, 1, .. , N-1} shifted indices ksm1: k-1 = {-1, 0, 1, .. , N-2} shifted indices ksp1: k+1 = {1, 2, .. , N-1, 0}

Type:

The Instanton is periodic

action()[source]

Action in au / fs, Hartree per femtosecond.

action_gradient()[source]

kin_grad corresponds to the gradient of S_0 in (Eq. 1 in [1], or first term in Eq. 6 in [2].) It boils down to the derivative of a sum of vector norms

d sum_k (||y_k - y_(k-1)||_2)² --- d y_k

The derivative of a norm of a vector difference is quite simple, but care has to be taken to recognize, that y_k appears two times in the sum. It appears in the first summand for k and in the second summand for k+1.

sum_k (||y_k - y_(k-1)||_2)²
  1. term 2. term

= (||y_k - y_(k-1)||_2)² + (||y_(k+1) - y_k||_2)² + ... and so on

The derivative of the first term is

2 * (y_k - y_(k-1))

and the derivative of the second term is

-2 * (y_(k+1) - y_k))

which is equal to

2 * (y_k - y_(k+1)) .

To summarize:

d sum_k(||y_k - y_(k-1)||_2)² --- d y_k

= 2 * (2 * y_k - y_(k-1) - y_(k+1)) .

action_hessian()[source]
as_xyz()[source]
property cart_coords
property cart_forces
property cart_hessian
property coords
property energy
property forces
classmethod from_instanton(other, **kwargs)[source]
classmethod from_ts(ts_geom, P, dr=0.4, delta_T=25, cart_hessian=None, **kwargs)[source]
get_additional_print()[source]
property gradient
property hessian
is_analytical_2d()[source]
property path_length
pysisyphus.irc.Instanton.T_crossover_from_eigval(eigval)[source]
pysisyphus.irc.Instanton.T_crossover_from_ts(ts_geom)[source]
pysisyphus.irc.Instanton.log_progress(val, key, i)[source]

17.1.1.10.11. pysisyphus.irc.LQA module

class pysisyphus.irc.LQA.LQA(geometry, N_euler=5000, **kwargs)[source]

Bases: IRC

step()[source]

17.1.1.10.12. pysisyphus.irc.ModeKill module

class pysisyphus.irc.ModeKill.ModeKill(geometry, kill_inds, nu_thresh=-5.0, do_hess=True, hessian_update=True, **kwargs)[source]

Bases: IRC

get_additional_print()[source]
postprocess()[source]
prepare(*args, **kwargs)[source]
step()[source]
update_mw_down_step()[source]

17.1.1.10.13. pysisyphus.irc.ParamPlot module

class pysisyphus.irc.ParamPlot.ParamPlot(coords_list, param1_inds, param2_inds)[source]

Bases: object

calc_angle(coords, ind1, ind2, ind3)[source]
calc_bond(coords, ind1, ind2)[source]
get_param(coords, param_inds)[source]
plot()[source]
save(path, prefix)[source]
save_coords(path, prefix)[source]
show()[source]

17.1.1.10.14. pysisyphus.irc.RK4 module

class pysisyphus.irc.RK4.RK4(geometry, step_length=0.1, max_cycles=125, downhill=False, forward=True, backward=True, root=0, hessian_init=None, displ='energy', displ_energy=0.001, displ_length=0.1, displ_third_h5=None, rms_grad_thresh=0.001, hard_rms_grad_thresh=None, energy_thresh=1e-06, imag_below=0.0, force_inflection=True, check_bonds=True, out_dir='.', prefix='', dump_fn='irc_data.h5', dump_every=5)[source]

Bases: IRC

get_k(mw_coords)[source]
step()[source]

17.1.1.10.15. pysisyphus.irc.initial_displ module

class pysisyphus.irc.initial_displ.ThirdDerivResult(coords_plus, energy_plus, H_plus, coords_minus, energy_minus, H_minus, G_vec, vec, ds)

Bases: tuple

G_vec

Alias for field number 6

H_minus

Alias for field number 5

H_plus

Alias for field number 2

coords_minus

Alias for field number 3

coords_plus

Alias for field number 0

ds

Alias for field number 8

energy_minus

Alias for field number 4

energy_plus

Alias for field number 1

vec

Alias for field number 7

pysisyphus.irc.initial_displ.cubic_displ(H, v0, w0, Gv, dE)[source]
According to Eq. (26) in [2] v1 does not depend on the sign of v0.

v1 = (F0 - 2v0^T F0 v0 I)⁻¹ x ([G0v0] - v0^T [G0v0] v0 I) v0

The first term is obviously independent of v0's sign. Using v0' = -v0 the second term becomes

([G0v0'] - v0'^T [G0v0'] v0' I) v0' (-[G0v0] - v0^T [G0v0'] v0 I) v0' (-[G0v0] + v0^T [G0v0] v0 I) v0' -(-[G0v0] + v0^T [G0v0] v0 I) v0 ([G0v0] - v0^T [G0v0] v0 I) v0

Strictly speaking the contraction [G0v0] should depend on the sign of v0. In the current implementation, the direction of v0 is not taken into account, but

get_curv_vec(H, Gv, v0, w0) == get_curv_vec(H, -Gv, -v0, w0) .

But somehow the Taylor expansion gives bogus results when called with -Gv and -v0...

pysisyphus.irc.initial_displ.cubic_displ_for_geom(geom, dE=-0.0005)[source]
pysisyphus.irc.initial_displ.cubic_displ_for_h5(h5_fn='third_deriv.h5', dE=-0.0005)[source]
pysisyphus.irc.initial_displ.get_curv_vec(H, Gv, v0, w0)[source]
pysisyphus.irc.initial_displ.taylor_closure(H, Gv, v0, v1, w0)[source]

Taylor expansion of energy to 3rd order.

dx(ds) = ds*v0 + ds**2*v1 / 2 dE(ds) = dx^T H dx / 2 + dx^T [Gv] dx / 6

H = Hessian Gv = 3rd derivative of energy along v0 v0 = Transition vector v1 = Curvature vector w0 = Eigenvalue belonging to v0

pysisyphus.irc.initial_displ.third_deriv_fd(geom, vec, ds=0.001)[source]

Third derivative of the energy in direction 'vec'.

17.1.1.10.16. Module contents

class pysisyphus.irc.DampedVelocityVerlet(geometry, v0=0.04, dt0=0.5, error_tol=0.003, max_cycles=150, **kwargs)[source]

Bases: IRC

damp_velocity(velocity)[source]
estimate_error(new_mw_coords)[source]
mw_grad_to_acc(mw_grad)[source]

Takes care of the units for the mass-weighted gradient. Converts units of a mass-weighted gradient [Hartree/(Bohr*amu)] to units of acceleration [sqrt(amu)*Bohr/fs²]. The 1e30 comes from converting second² to femto second².

prepare(direction)[source]
step()[source]
class pysisyphus.irc.Euler(geometry, step_length=0.01, **kwargs)[source]

Bases: IRC

step()[source]
class pysisyphus.irc.EulerPC(*args, hessian_recalc=None, hessian_update='bofill', hessian_init='calc', max_pred_steps=500, loose_cycles=3, dump_dwi=False, scipy_method=None, corr_func='mbs', **kwargs)[source]

Bases: IRC

corrector_step(init_mw_coords, step_length, dwi)[source]
get_integration_length_func(init_mw_coords)[source]
prepare(*args, **kwargs)[source]
scipy_corrector_step(init_mw_coords, step_length, dwi)[source]

Solve IRC equation dx/ds = -g/|g| on DWI PES in mass-weighted coordinates. Integration done until self.step_length in unweighted coordinates is achieved.

step()[source]
class pysisyphus.irc.GonzalezSchlegel(geometry, max_micro_cycles=20, micro_step_thresh=0.001, hessian_recalc=None, line_search=False, **kwargs)[source]

Bases: IRC

micro_step(counter)[source]

Constrained optimization on a hypersphere.

perp_component(vec, perp_to)[source]
postprocess()[source]
step()[source]
class pysisyphus.irc.IMKMod(geometry, corr_first=True, corr_first_size=0.5, corr_first_energy=True, corr_bisec_size=0.0025, corr_bisec_energy=True, **kwargs)[source]

Bases: IRC

fit_grad_and_energies(energy0, grad0, energy1, direction)[source]
fit_parabola(x, y)[source]
step()[source]
class pysisyphus.irc.LQA(geometry, N_euler=5000, **kwargs)[source]

Bases: IRC

step()[source]
class pysisyphus.irc.ModeKill(geometry, kill_inds, nu_thresh=-5.0, do_hess=True, hessian_update=True, **kwargs)[source]

Bases: IRC

get_additional_print()[source]
postprocess()[source]
prepare(*args, **kwargs)[source]
step()[source]
update_mw_down_step()[source]
class pysisyphus.irc.RK4(geometry, step_length=0.1, max_cycles=125, downhill=False, forward=True, backward=True, root=0, hessian_init=None, displ='energy', displ_energy=0.001, displ_length=0.1, displ_third_h5=None, rms_grad_thresh=0.001, hard_rms_grad_thresh=None, energy_thresh=1e-06, imag_below=0.0, force_inflection=True, check_bonds=True, out_dir='.', prefix='', dump_fn='irc_data.h5', dump_every=5)[source]

Bases: IRC

get_k(mw_coords)[source]
step()[source]