GULP¶
GULP, the General Utility Lattice Program, is a forcefield code.
Instructions¶
Make sure that the following variables are set correctly:
$GULP_LIB: path to the folder containing the potential files
$ASE_GULP_COMMAND="/path/to/gulp < PREFIX.gin > PREFIX.got"
The latter defaults to gulp < PREFIX.gin > PREFIX.got, assuming that
the gulp command is in $PATH.
GULP uses several variables to write the input file. The most important ones are:
keywords: string of space-separated keywords, e.g.:'angle compare comp conp'
options: list of strings with options, e.g.:['xtol opt 5.0']
library: filename of the potential library file (e.g.'reaxff.lib'). Recall that the search path for this library is given by$GULP_LIB.
If the potential uses different atom types, one must define an instance of conditions
with the atom object to label them and add rules to rename atoms. The only rule suported
right now is min_distance_rule which selects atoms of one type close to atoms of another type and
renames them. For example:
c = Conditions(atoms)
c.min_distance_rule('O', 'H', ifcloselabel1='O2', ifcloselabel2='H',
elselabel1='O1')
calc = GULP(conditions=c)
This will assign, for each H in the system, a corresponding O atom as O2, selecting those O which are the closest to H. The rest of oxygens will be labeled O1.
Finally, if the potential file requires the use of shells, the
variable shel must be used:
GULP(..., shel=['O'])
Example¶
Here is an example of how to use the GULP calculator.
# flake8: noqa
import numpy as np
from ase import Atoms
from ase.calculators.gulp import GULP, Conditions
cluster = Atoms(symbols='O4SiOSiO2SiO2SiO2SiOSiO2SiO3SiO3H8',
pbc=np.array([False, False, False], dtype=bool),
cell=np.array(
[[ 0., 0., 0.],
[ 0., 0., 0.],
[ 0., 0., 0.]]),
positions=np.array(
[[-1.444348, -0.43209 , -2.054785],
[-0.236947, 2.98731 , 1.200025],
[ 3.060238, -1.05911 , 0.579909],
[ 2.958277, -3.289076, 2.027579],
[-0.522747, 0.847624, -2.47521 ],
[-2.830486, -2.7236 , -2.020633],
[-0.764328, -1.251141, 1.402431],
[ 3.334801, 0.041643, -4.168601],
[-1.35204 , -2.009562, 0.075892],
[-1.454655, -1.985635, -1.554533],
[ 0.85504 , 0.298129, -3.159972],
[ 1.75833 , 1.256026, 0.690171],
[ 2.376446, -0.239522, -2.881245],
[ 1.806515, -4.484208, -2.686456],
[-0.144193, -2.74503 , -2.177778],
[ 0.167583, 1.582976, 0.47998 ],
[-1.30716 , 1.796853, -3.542121],
[ 1.441364, -3.072993, -1.958788],
[-1.694171, -1.558913, 2.704219],
[ 4.417516, 1.263796, 0.563573],
[ 3.066366, 0.49743 , 0.071898],
[-0.704497, 0.351869, 1.102318],
[ 2.958884, 0.51505 , -1.556651],
[ 1.73983 , -3.161794, -0.356577],
[ 2.131519, -2.336982, 0.996026],
[ 0.752313, -1.788039, 1.687183],
[-0.142347, 1.685301, -1.12086 ],
[ 2.32407 , -1.845905, -2.588202],
[-2.571557, -1.937877, 2.604727],
[ 2.556369, -4.551103, -3.2836 ],
[ 3.032586, 0.591698, -4.896276],
[-1.67818 , 2.640745, -3.27092 ],
[ 5.145483, 0.775188, 0.95687 ],
[-2.81059 , -3.4492 , -2.650319],
[ 2.558023, -3.594544, 2.845928],
[ 0.400993, 3.469148, 1.733289]]))
c = Conditions(cluster)
c.min_distance_rule('O', 'H', ifcloselabel1='O2',
ifcloselabel2='H', elselabel1='O1')
calc = GULP(keywords='conp', shel=['O1', 'O2'], conditions=c)
# Single point calculation
cluster.calc = calc
print(cluster.get_potential_energy())
# Optimization using the internal optimizer of GULP
calc.set(keywords='conp opti')
opt = calc.get_optimizer(cluster)
opt.run(fmax=0.05)
print(cluster.get_potential_energy())
The script performs a single-point calculation and then an optimization
of the H8Si8O20 double ring
structure within ASE using the ffsioh potential. The
method get_optimizer() returns an object which works like an ASE optimizer,
but actually triggers an optimization using GULP’s internal optimizer.
This works on GULP calculators that have the 'opti' keyword.