# Structure and Reactivity¶

To study structure and reactivity common tasks are **geometry optimization**, **transition state search** and a **frequency run**. A frequency run can be done to obtain a good initial Hessian for a geometry optimization, or to check whether a proper minimum or transition state has been found.

The basis for locating stationary points on the PES are the nuclear gradients. (The user may assume that the program does that automatically when needed.)

## Nuclear energy gradients¶

If you want the program to calculate the gradient of the energy w. r. t. nuclear displacements you should add the `Gradients`

keyword. If the `SelectedAtoms`

key is present, only the gradient of the selected atoms will be evaluated, and the rest will be set to zero.

Details about the theory of the analytical gradients within the BAND program can be found in the work of *Kadantsev* et al. [19].

`Gradients (block-type)`

Key to control the geometry optimization:

Gradients {LatticeGradients [false | true]} {LatticeStepSize step} {UseRelativeLatticeStep [true | false]} End

`LatticeGradients`

- (
**Default: false**) This key handles whether the lattice gradients are evaluated numerically (false) or analytically (true). `latticeStepSize`

- (
**Default: 0.012**) With this keyword the step size for the numerical differentiation can be changed. When`UseRelativeLatticeStep`

is set to true the step will be relative to the size of the displacement. (Unit is Bohr)

The user rarely needs to specify the `Gradients`

keyblock as it is automatically added for geometry optimizations and frequency runs.

This option honors the `SelectedAtoms`

key, in which case only the forces will be calculated for the selected atoms.

## Lattice gradients¶

The step size for numerical differentiation is controlled with Gradients%latticeStepSize, and Gradients%UseRelativeLatticeStep.

`LatticeGradients (block-type)`

Key to handle special options for the calculation of analytical and numerical lattice gradients.

LatticeGradients {Analytical [true | false]} {AnalyticalPulay [true | false]} {AnalyticalKinetic [true | false]} {AnalyticalXC [true | false]} {AnalyticalElectrostatic [true | false]} End

`Analytical`

- (
**Default: false**) If set to true, the following options`AnalyticalElectrostatic`

,`AnalyticalKinetic`

,`AnalyticalPulay`

, and`AnalyticalXC`

will be set to true automatically. `AnalyticalElectrostatic`

- (
**Default: true**) If true the electrostatic energy contribution to the lattice gradients is calculated analytically, else it is calculated numerically. `AnalyticalKinetic`

- (
**Default: true**) If true the kinetic energy contribution to the lattice gradients is calculated analytically, else it is calculated numerically. `AnalyticalPulay`

- (
**Default: true**) If true the Pulay contribution to the lattice gradients is calculated analytically, else it is calculated numerically. `AnalyticalXC`

- (
**Default: true**) If true the exchange-correlation energy contribution to the lattice gradients is calculated analytically, else it is calculated numerically.

In case you want to calculate forces acting on lattice vectors without performing an geometry optimization, you should use `LatticeGradients%Analytical true`

!

## Geometry optimization (GeoOpt)¶

The geometry can be optimized by adding the `GeoOpt`

keyblock. Within the block you can specify the maximum number of cycles and the convergence criterion. The optimizer works with Cartesian coordinates and a Quasi Newton stepper. The initial Hessian is by default the unit matrix, but can also be loaded from a previous geometry optimization or frequency run, also from faster methods as DFTB and UFF. Various constraints are possible such as bond distances and angles. The GDIIS convergence accelerator is available but is disabled by default. A geometry optimization can be restarted from a previous result file. The geometry can be optimized non-relativistically, with scalar relativistic ZORA, or with spin-orbit coupling.

See also

Troubleshooting: Geometry does not converge

`GeoOpt (block-type)`

Key to control the automatic geometry optimization:

GeoOpt {OptimizeLattice [false | true]} {Iterations n} {Converge {Grad=tolgrad} {E=tolE} {Step=tolStep}} {TrustRadius radius} {UseVariableTrustRadius var} {InitialHessian [UnitMatrix | filename | Swart]} {RestartSCF [true | false]} {RestartFit [true | false]} {GDIIS [true | false]} {StaticCosmoSurface [true | false]} {TSMode mode} End

`OptimizeLattice`

- (
**Default: false**) Whether or not to optimize the lattice vectors. The lattice gradients are obtained by numerical or analytical differentiation. `Iterations`

- (
**Default: 50**) maximum number of cycles. `Converge {Grad=tolgrad} {E=tolE} {Step=tolStep}`

- various criteria to determine the convergence can be specified on line. Here
*tolGrad*is the maximum gradient allowed (**Default: 1.0e-2**Hartree per Angstrom),*tolE*is the maximum energy change allowed (**Default: 1.0e-3**Hartree),*tolStep*is the maximum step size (**Default: 3e-2**Bohr) `TrustRadius`

- (
**Default: 0.2**Bohr) The step sizes taken by the optimizer will be limited to this value. If the proposed step is larger, it will be scaled back. `UseVariableTrustRadius`

- (
**Default: false**) Automatic adjustment of the trust radius based on the observed energy changes. `InitialHessian`

- (
**Default: UnitMatrix**) The initial Hessian is read from a “RUNKF” file*filename*of a previous frequency calculation. Band can also read the Hessian from DFTB or UFF calculations (either a frequency run or a geometry optimization). Another option,*Swart*, is to use the empirical force field derived Swart Hessian. `RestartSCF`

- (
**Default: true**) During the optimization try to restart the SCF from the previous geometry step. `RestartFit`

- (
**Default: false**) During the optimization try to restart the SCF from fit coefficients of the previous geometry step. Can not be used simultaneously with RestartSCF. `GDIIS`

- (
**Default: false**) Use the GDIIS convergence accelerator. `StaticCosmoSurface`

- (
**Default: false**) Keep the Cosmo surface (if any) constant. `TSMode`

- Hessian normal mode index to follow during transition state search. The specified mode is selected on the first iteration and is followed regardless of its index on subsequent iterations. Modes with zero frequencies are not counted.

Note that a (failed) geometry optimization can be continued with the Restart key by specifying as option GeometryOptimization.

## Convergence¶

A geometry optimization is considered *converged* when all the following criteria are met (this information is also printed to the logfile):

- the difference between the bond energy at the current geometry and at the previous geometry step is smaller than \(\text{TolE}\)
- the maximum Cartesian nuclear gradient is smaller than \(\text{TolG}\)
- the root mean square (RMS) of the Cartesian nuclear gradients is smaller than \(\frac{2}{3} \text{TolG}\)
- the maximum Cartesian step is smaller than \(\text{TolR}\)
- the root mean square (RMS) of the Cartesian steps is smaller than \(\frac{2}{3} \text{TolR}\)

**Note**: if the maximum and RMS gradients are 10 times smaller than the convergence criterion, then criteria 4 and 5 are ignored

The convergence criteria can be changed via the `Converge`

key:

```
GEOMETRY
Converge {E=TolE} {Grad=TolG} {Rad=TolR} {Angle=TolA}
END
```

`TolE`

- The criterion for changes in the energy, in Hartree.
**Default: 1e-3**. `TolG`

- The criterion for the nuclear gradients, in Hartree/Angstrom.
**Default: 1e-3**. `TolR`

- Refers to changes in the Cartesian coordinates or bond lengths, depending on in what coordinates you optimize, in Angstrom.
**Default: 1e-2**. `TolA`

- Refers to changes in bond- and dihedral angles, in degrees. This is only meaningful if optimization takes place in Z-matrix coordinates.
**Default: 0.5 degree**.

If only a numerical value is supplied as argument for `Converge`

, rather than a specification by name, it is considered to apply to the gradients (only):

```
GEOMETRY
Converge 1E-4
End
```

is equivalent to:

```
GEOMETRY
Converge grad=1E-4
End
```

## Numerical frequencies (Hessian)¶

By specifying the `RunType`

*Frequencies* you can calculate the Hessian which in turn yields the harmonic frequencies and vibrational modes (*q → 0 limit* of a phonons calculation). The Hessian is an important product that can be used in a geometry optimization, and in particular to start a transition state search. The second derivative of the energy (i.e. the Hessian) is obtained by numerically differentiating the analytical gradients. When your unit cell has N atoms in the unit cell and no symmetry (other than translational symmetry) the number of displacements is 2x3xN. For large systems this can be very time consuming, and you could consider to calculate only a partial Hessian. If you have a very symmetric unit cell, but are interested in asymmetric modes, set the subkeys `useA1Displacements`

and `doA1Projection`

to false. **Note**: the lattice vectors are assumed to be fixed.

A frequency run is invoked with

```
RunType
Frequencies
End
```

The rest is controlled by the `Frequencies`

key.

`Frequencies (block-type)`

Key to control the numerical frequency run:

Frequencies {Step size} {useA1Displacements [true | false]} {doA1Projection [true | false]} {StaticCosmoSurface [true | false]} End

`Step`

- (
**Default: 0.01**) The step size (in Å) for the numerical differentiation. `useA1Displacements`

- (
**Default: true**) Determine only the total symmetric modes. When disabled, NOSYM calculations will be performed. `doA1Projection`

- (
**Default: true**) Symmetrize the Hessian. For most users the only sensible setting is to make it equal to the value of useA1Displacements. `StaticCosmoSurface`

- (
**Default: true**) Keep Cosmo surface (if any) constant. Although this is an approximation it is numerically much more stable. If you set it to*false*there is the risk that the number of Cosmo points is not constant for all finite displacements.

This option honors the

`SelectedAtoms`

key, in which case only the Hessian will be calculated for the selected atoms.

## Transition state search¶

A transition state search is technically a geometry optimization. It is a search in the 3xN space for a point with vanishing gradients. The only difference is that it will try to go uphill in the direction of the lowest vibrational mode. The normal procedure is to guess a point that is as close as possible to the transition state. In that point you do a frequency run, which should hopefully produce a lowest vibrational mode in the direction of the transition state. You then start the Transition state search from that point, using the Hessian of the frequency run. At the end you can do again a frequency run to check that there is indeed only one single negative vibrational mode.

The transition state search is invoked like this

```
RunType
TS
End
```

The rest is controlled by the `GeoOpt`

key as in a normal geometry optimization.

## Partial Hessian and (pre-)optimizations¶

If you consider the interaction of a molecule with a surface it may be initially convenient to freeze some of the lower layers in the slab, making calculations much cheaper. This can be done for the evaluation of the Hessian (see SelectedAtoms) and during a geometry optimization or transition state search (see Constraints).

Caution

When using a partial Hessian to start a geometry optimization note that the atoms that you select during the frequency run, are the ones that should not be frozen in the geometry optimization. In a way the selection needs to be inverted.

## Constrained optimization¶

During a geometry optimization certain constraints can be enforced with the `Constraints`

key. The constraints refer to the coordinates of the unit cell as specified on input.

`Constraints (block-type)`

Key to control constraints during a geometry optimization:

Constraints {Atom aI {x y z} {Coord aI coordI {coordIvalue} {Dist aI aJ distance} {Angle aI aJ aK angle} {Dihed aI aJ aK aL angle} End

`Atom`

- Freeze the position of atom
*a1*to the current coordinate, or if specified to (*x y z*). This constraint will make the calculation cheaper, because the gradients of frozen atoms need not be calculated. `Dist`

- Constrain the distance between atoms
*a1*and*a2*to*distance*. `Angle`

- Constrain the angle between atoms
*a1*,*a2*, and*a3*to*angle*. `Dihed`

- Constrain the dihedral angle between atoms
*a1*,*a2*,*a3*, and*a4*to*angle*.

## Selected atoms¶

When doing a frequency or geometry optimization run you can restrict it to a limited number of selected nuclei. The selection may not break the symmetry

```
SelectedAtoms a1 a2 ... an
```

The list of atoms with indices ‘a1 a2 ... an’ are considered selected.

**Note**: Currently this option only has effect in a numerical frequency run and with the Gradients keyword.

## Run Types (GeoOpt, TS, Frequencies and Phonons)¶

`RunType (block-type)`

Key to invoke the following calculation types:

RunType GeoOpt End

RunType TS End

RunType Frequencies End

RunType Phonons End