Geometry Optimisation

DL-FIND is a powerful and flexible geometry optimisation library and is the standard optimiser module in ChemShell. The code is under active development but the features documented here are those that are available in the version distributed with ChemShell.

DL-FIND supports a wide variety of optimisation methods. Standard local minimisation techniques are supported (steepest descent, conjugate gradient, Newton-Raphson, BFGS) as well as damped dynamics and limited memory BFGS (L-BFGS). Transition state methods include nudged elastic band (NEB) for finding reaction paths and P-RFO and the dimer method for refining transition states. Conical intersection optimisations are supported using penalty function, gradient projection and Lagrange-Newton algorithms. Parallel optimisation methods include a genetic algorithm, stochastic search and parallel NEB.

DL-FIND accepts input structures in Cartesian coordinates. These may be converted into a number of internal coordinate systems for the optimisation, including delocalised internal coordinates (DLCs), total connection, and hybrid delocalised internal coordinates (HDLCs) for large systems. Constrained coordinates are supported as well as a general mechanism for imposing restraints.

A detailed description of DL-FIND’s features can be found in:

Johannes Kästner, Joanne M. Carr, Thomas W. Keal, Walter Thiel, Adrian Wander and Paul Sherwood, J. Phys. Chem. A, 2009, 113 (43), 11856-11865. DOI: 10.1021/jp9028968

Please cite this paper if you use DL-FIND in published work.

Running an optimisation

class Opt

In ChemShell, DL-FIND optimisations can be run using the Opt() driver, which expects a Theory object, e.g.:

my_theory = NWChem(...)
my_job = Opt(theory=my_theory)

The geometry optimisation calculation is run with:

my_job.run()

Then the results will be stored in my_job.result (see: <class Result>).

Methods

Opt.run(dryrun=False)

Parameters: dryrun (bool, optional) – If the function is called with dryrun=True, the job will only generate input scripts without executing real calculatons

General options

DL-FIND optimisations take the following arguments:

theory: Choose the theory level by specifying a Theory instance

frag2: Initial second configuration of the system. End point of the initial path in NEB calculations, direction of the dimer in dimer calculations. Can also be a list of fragments, in which case it defines points along the initial NEB path (the desired endpoint should be given last).

maxcycle: (default: 100) Maximum number of optimisation cycles

maxene: (default: 10000) Maximum number of energy evaluations

Coordination system options

coordinates

Specifies which coordinate system to use

Allowed values:

'cartesian': (default) Cartesian coordinates
'mass': Mass-weighted Cartesian coordinates
'dlc': Delocalised internal coordinates (DLC). A redundant set of bonds, angles, torsions, and sometimes inversions is created. A non-redundant combination of them is found by diagonalising the spectroscopic G-matrix. The optimisation is performed in this non-redundant set. `
'tc': Same as 'dlc', but the redundant set consists only of bond - all atoms in the system are connected (total connection scheme)
'hdlc': Hybrid delocalised internal coordinates. The system is partitioned into fragments. Delocalised coordinates (as in DLC) are used within each fragment. The fragments are coupled via Cartesian coordinates. This version is recommended for large systems, as the coordinate transformation scales linearly with the number of atoms.
'hdlc-tc': Hybrid delocalised internal coordinates based on the total connection scheme.

NB: that the 'dlc' option is not equivalent to specifying hdlc with a single residue. The 'dlc' case gives pure internal delocalised coordinates (6 external degrees of freedom are dropped), whereas with 'hdlc' Cartesian coordinates are added to the primitive list even if there is only one residue, and the external degrees of freedom are retained.

residues

Defines residues (atom groups) for use with HDLC coordinates.

Accepts a dictionary of the form:

{ 'res1' : [atom1, atom2, atom3, ...],
  'res2' : [atom4, atom5, atom6, ...] }

NB: For biomolecular Fragments this dictionary can be generated from PDB data using the getHDLCResidues() method on the fragment.

constraints

Define constraints for use with either HDLC or Cartesian coordinates.

For HDLC coordinates, internal constraints can be defined using a list of the form:

[ ('bond', atom1, atom2),
  ('angle', atom1, atom2, atom3),
  ('torsion', atom1, atom2, atom3, atom4),
 ... ]

For Cartesian coordinates, constraints can be defined using a list of the form:

[ ('x', atom1), ('yz', atom2), ... ]

where the first constraint freezes the X coordinate of atom1, the second freezes the Y and Z coordinates of atom2, etc.

Optimisation algorithm options

algorithm

Specifies the optimisation algorithm

Allowed values:

'lbfgs': (default) Limited-memory Broyden–Fletcher–Goldfarb–Shanno (LBFGS) optimisation. Recommended for minimum search, NEB, and the dimer method. The time and memory requirements spent for determining the search direction and step length scale linearly with the system size. The number of steps kept in memory is specified by lbfgs_mem.
'prfo': The partitioned rational function optimisation method. Recommended for transition state search (unless done by the dimer or NEB methods). It requires the calculation of the Hessian. The keyword update_method can be used to specify how the Hessian should be updated.
'cg': Conjugate gradient method following Polak–Ribiere.
'sd': Steepest descent.
'dyn': Damped molecular dynamics. Four additional parameters can be specified. timestep is the time step in atomic units. fric0 is the start friction, fricfac is the factor to reduce the friction each time the energy decreases, and fricp is the friction to apply if the energy increases. The frictions are defined so that 0 corresponds to free (undamped) dynamics, and 1 corresponds to steepest descent. This variable friction facilitates convergence to an energy minimum.
'nr': Basic Newton-Raphson/quasi-Newton optimiser.
'ln': Lagrange-Newton constrained optimiser (for conical intersection searches only).

trustradius

Specifies the trust radius approach

Allowed values:

'constant': (default) Constant trust radius (equal to maxstep)
'energy': Trust radius based on the energy change
'gradient': Trust radius based on the projection of the gradient on the current step. Can be used when no energy expression is available (e.g. NEB).

maxstep: (Default: 0.5) Maximum change of one coordinate component in one optimisation step.

tolerance: (Default: 0.00045) Convergence criterion for the maximum gradient component. The criterion for the maximum step component is 4*tolerance, for the RMS step it is 8/3*tolerance, for the RMS gradient it is 2/3*tolerance. All criteria are applied in the coordinate system specified by coordinates.

tolerance_e: (Default: tolerance/450) Convergence criterion for the energy change (in Hartree).

Nudged elastic band options

neb

Specifies if and how the nudged elastic band (NEB) method is to be used

Allowed values:

'no': (default) No NEB calculation will be performed. neb is automatically set to 'no' if no frag2 is specifed.
'free': The endpoints of the NEB path are minimised
'frozen': The endpoints of the NEB path are frozen
'perpendicular': The endpoints of the NEB path minimised perpendicular to the NEB path. Their distance to the nearest other point is kept fixed.

nimages: (Default: 1) Specifies the number of images in a NEB calculation.

nebk: (Default: 0.01) Force constant for NEB calculations.

Dimer method options

dimer

Specifies if the dimer method for transition state (TS) search is to be used.

Allowed values:

False: (default) Do not use the dimer method
True: Use the dimer method for TS search

delta: (Default: 0.01) Dimer separation (in internal coordinates)