Geometry optimization
CASTEP supports a number of schemes for geometry optimization: BFGS, TPSD, and damped molecular dynamics.
BFGS geometry optimization
The main advantage of the BFGS minimizer (Pfrommer et al., 1997) is the ability to perform cell optimization, including optimization at fixed external stress. The BFGS scheme uses a starting Hessian that the optimization updates recursively. The CASTEP implementation involves a Hessian in the mixed space of internal and cell degrees of freedom, allowing optimization of both lattice parameters and atomic coordinates.
The LBFGS algorithm (low-memory BFGS) is expected to perform better for large systems. This approach uses only a limited number of inverse Hessian updates to construct a new Hessian. This scheme scales linearly in memory, while BFGS scales as a square of the system size (Aarons, 2011). The LBFGS algorithm also implements a universal sparse preconditioner that accelerates geometry optimization by a factor of two for small systems and a factor of ten for very large ones (Packwood et al., 2016).
When you require cell optimization, use a finite basis set correction term. Calculating the term is not too costly (10-30% of the self-consistent electronic minimization at the first iteration) relative to the advantages that it provides. Cell optimization runs may be problematic if the finite basis set correction is not used or if the energy cutoff is so low that the correction is not accurate. In these circumstances, the optimization run may stop with a message that states that the energy has converged, but that the stress is still non-zero. The minimizer attempts to find the energy minimum rather than the stress zero-point, since the former is more meaningful under the circumstances.
Another potential pitfall is the use of the Fixed Basis Size Cell optimization setting when the starting geometry is very different from the final one. The finite basis set correction depends on the cell size and shape, although the minimizer disregards this dependence. When using Fixed Basis Size modifies the cell geometry, the effective cutoff energy changes but the number of plane waves remains fixed. If this change takes the dEtot/d(ln Ecut) function far away from the point used to evaluate the finite basis set correction, the results are not accurate. Therefore, compare the starting and final geometries and perform a completely new run starting from the final configuration if the difference between the two structures is large.
Constraints
CASTEP can perform geometry optimization with constraints applied. The simplest type of constraint is to fix the atom positions. In CASTEP, this means fixing the fractional coordinates of the atoms.
You do not need to manually fix atoms on special positions as the symmetry treatment in CASTEP ensures the correct behavior of such atoms.
It is also possible to fix lattice parameters. You do not need to fix lattice parameters, except for parameters not fixed by symmetry. This type of constraint is often useful in phase stability studies when, for example, you want to fix the angles to 90° (Winkler and Milman, 1997).
It is also possible to impose more general linear constraints on atomic coordinates. A matrix that transforms the Cartesian coordinates of all atoms to the subspace of unconstrained coordinates specifies these constraints. You can use this facility to fix, for example, the z-coordinates of atoms in slab calculations of surface processes.
Non-linear constraints
Non-linear constraints refer to constraints on interatomic distances (bonds), angles, and torsions. You can impose such constraints using the delocalized internals optimizer.
Estimated compressibility
You can facilitate cell optimization convergence by making a reasonable initial estimate of the bulk modulus for the material. The initial changes to lattice vectors are calculated on the basis of this value, based on the calculated stress tensor. When the estimated bulk modulus is too high, the steps become smaller than is ideal for the minimizer. Therefore, it is conceivably possible that you could reduce the number of geometry optimization steps by reducing this value. Conversely, too small a bulk modulus for a hard material could result in an oscillating behavior, where the first step is too big. Thus, prior knowledge of the hardness of the material could help to reduce the number of steps required for calculations. The table below illustrates the effect of a judicious choice of the initial bulk modulus value. The times are relative to the same calculation using an initial value of 500 GPa.
| Compound | Estimated bulk modulus (GPa) | Number of steps | Relative time |
|---|---|---|---|
| Cubane | 5.000 | 15 | 0.98 |
| 50.00 | 15 | 0.87 | |
| 250.0 | 15 | 0.82 | |
| 500.0 | 20 | 1.00 | |
| 1000.0 | 24 | 1.29 | |
| Glycine | 5.000 | 46 | 1.29 |
| 50.00 | 42 | 1.03 | |
| 250.0 | 39 | 0.92 | |
| 500.0 | 37 | 1.00 | |
| 1000.0 | 54 | 1.24 | |
| YBa2Cu3O7 | 5.000 | 12 | 0.81 |
| 50.00 | 13 | 0.74 | |
| 250.0 | 12 | 0.81 | |
| 500.0 | 17 | 1.00 | |
| 1000.0 | 13 | 0.81 |
TPSD geometry optimization
The BFGS algorithm as implemented in CASTEP uses the assumption that the potential is quadratic to determine an optimal step length per iteration. When this is not the case, the BFGS algorithm does not perform well and may take a long time to optimize the structure. A more robust Barzilai-Borwein (1988) Two-Point Steepest Descent (TPSD) algorithm allows you to solve a larger range of problems. It is specifically more efficient when you request cell optimization with custom constraints. This scenario describes, for example, optimization of solid-solid interfaces strained to a substrate: that is, with fixed lateral cell parameters and all cell angles, only allowing relaxation of the cell dimension in the direction normal to the interface.
The TPSD algorithm requires only the gradient and position of the system at the current and previous iteration, meaning that the calculation consumes a relatively small amount of memory. The difference of each value between the iterations generates a step size for use in the direction of the most negative gradient from the current iteration. TPSD optimization converges reliably in a number of situations where the BFGS algorithm takes too long to converge.
Damped molecular dynamics
Damped molecular dynamics presents an alternative method for geometry minimization that involves only internal coordinates (with fixed cell parameters). The method uses the critical damping regime as a way of dealing with the ground state. The regime can be implemented either by using one damping coefficient for all degrees of freedom (coupled modes) or by using different coefficients for different degrees of freedom (independent modes). The latter approach allows you to freeze all modes, both fast and slow, simultaneously. Alternatively, you can perform steepest descents damping with a fixed coefficient. However, this is a less efficient approach. Indeed, it is not, strictly speaking, a molecular dynamics technique, since it solves first-order equations of motion and not second-order equations.
Both independent modes and coupled modes damped molecular dynamics runs can be performed with a bigger time step than undamped molecular dynamics simulations. When the system gets close to equilibrium, the time step can be increased even further, since the fast modes freeze out before the slow ones. CASTEP can automatically adjust the time step, leading to increased efficiency of the algorithm. It is also recommended that you recalculate the damping coefficients periodically during a damped molecular dynamics run. A full description of the method implemented in CASTEP can be found elsewhere (Probert, 2003).