CASTEP [1] W is a software package which uses density functional theory to provide a good atomic-level description of all manner of materials and molecules. CASTEP can give information about total energies, forces and stresses on an atomic system, as well as calculating optimum geometries, band structures, optical spectra, phonon spectra and much more. It can also perform molecular dynamics simulations.
The CASTEP calculator interface class offers intuitive access to all CASTEP settings and most results. All CASTEP specific settings are accessible via attribute access (i.e. calc.param.keyword = ... or calc.cell.keyword = ...)
Set the environment variables appropriately for your system.
>>> export CASTEP_COMMAND=' ... '
>>> export CASTEP_PP_PATH=' ... '
>>> export PSPOT_DIR=' ... '
The default initialization command for the CASTEP calculator is
To do a minimal run one only needs to set atoms, this will use all default settings of CASTEP, meaning LDA, singlepoint, etc..
With a generated castep_keywords.py in place all options are accessible by inspection, i.e. tab-completion. This works best when using ipython. All options can be accessed via calc.param.<TAB> or calc.cell.<TAB> and documentation is printed with calc.param.<keyword> ? or calc.cell.<keyword> ?. All options can also be set directly using calc.keyword = ... or calc.KEYWORD = ... or even calc.KeYwOrD or directly as named arguments in the call to the constructor (e.g. Castep(task='GeometryOptimization')).
All options that go into the .param file are held in an CastepParam instance, while all options that go into the .cell file and don’t belong to the atoms object are held in an CastepCell instance. Each instance can be created individually and can be added to calculators by attribute assignment, i.e. calc.param = param or calc.cell = cell.
All internal variables of the calculator start with an underscore (_). All cell attributes that clearly belong into the atoms object are blocked. Setting calc.atoms_attribute (e.g. = positions) is sent directly to the atoms object.
Keyword | Description |
---|---|
directory | The relative path where all input and output files will be placed. If this does not exist, it will be created. Existing directories will be moved to directory-TIMESTAMP unless self._rename_existing_dir is set to false. |
label | The prefix of .param, .cell, .castep, etc. files. |
Internal Setting | Description |
---|---|
_castep_command | (=castep): the actual shell command used to call CASTEP. |
_check_checkfile | (=True): this makes write_param() only write a continue or reuse statement if the addressed .check or .castep_bin file exists in the directory. |
_copy_pspots | (=False): if set to True the calculator will actually copy the needed pseudo-potential (*.usp) file, usually it will only create symlinks. |
_export_settings | (=True): if this is set to True, all calculator internal settings shown here will be included in the .param in a comment line (#) and can be read again by merge_param. merge_param can be forced to ignore this directive using the optional argument ignore_internal_keys=True. |
_force_write | (=True): this controls wether the *cell and *param will be overwritten. |
_prepare_input_only | (=False): If set to True, the calculator will create *cell und *param file but not start the calculation itself. If this is used to prepare jobs locally and run on a remote cluster it is recommended to set _copy_pspots = True. |
_castep_pp_path | (='.') : the place where the calculator will look for pseudo-potential files. |
_rename_existing_dir | (=True) : when using a new instance of the calculator, this will move directories out of the way that would be overwritten otherwise, appending a date string. |
_set_atoms | (=False) : setting this to True will overwrite any atoms object previously attached to the calculator when reading a .castep file. By de- fault, the read() function will only create a new atoms object if none has been attached and other- wise try to assign forces etc. based on the atom’s positions. _set_atoms=True could be necessary if one uses CASTEP’s internal geometry optimization (calc.param.task='GeometryOptimization') because then the positions get out of sync. Warning: this option is generally not recommended unless one knows one really needs it. There should never be any need, if CASTEP is used as a single-point calculator. |
_track_output | (=False) : if set to true, the interface will append a number to the label on all input and output files, where n is the number of calls to this instance. Warning: this setting may con- sume a lot more disk space because of the additio- nal *check files. |
_try_reuse | (=_track_output) : when setting this, the in- terface will try to fetch the reuse file from the previous run even if _track_output is True. By de- fault it is equal to _track_output, but may be overridden. Since this behavior may not always be desirable for single-point calculations. Regular reuse for e.g. a geometry-optimization can be achieved by setting calc.param.reuse = True. |
[1] | S. J. Clark, M. D. Segall, C. J. Pickard, P. J. Hasnip, M. J. Probert, K. Refson, M. C. Payne Zeitschrift für Kristallographie 220(5-6) pp.567- 570 (2005) PDF. |
"""This simple demo calculates the total energy of CO molecules
using once LDA and once PBE as xc-functional. Obviously
some parts in this scripts are longer than necessary, but are shown
to demonstrate some more features."""
import ase
import ase.calculators.castep
import ase.io.castep
calc = ase.calculators.castep.Castep()
directory = 'CASTEP_ASE_DEMO'
# include interface settings in .param file
calc._export_settings = True
# reuse the same directory
calc._directory = directory
calc._rename_existing_dir = False
calc._label = 'CO_LDA'
# necessary for tasks with changing positions
# such as GeometryOptimization or MolecularDynamics
calc._set_atoms = True
# Param settings
calc.param.xc_functional = 'LDA'
calc.param.cut_off_energy = 400
# Prevent CASTEP from writing *wvfn* files
calc.param.num_dump_cycles = 0
# Cell settings
calc.cell.kpoint_mp_grid = '1 1 1'
calc.cell.fix_com = False
calc.cell.fix_all_cell = True
# Set and clear and reset settings (just for shows)
calc.param.task = 'SinglePoint'
# Reset to CASTEP default
calc.param.task.clear()
# all of the following are identical
calc.param.task = 'GeometryOptimization'
calc.task = 'GeometryOptimization'
calc.TASK = 'GeometryOptimization'
calc.Task = 'GeometryOptimization'
# Prepare atoms
mol = ase.atoms.Atoms('CO', [[0, 0, 0], [0, 0, 1.2]], cell=[10, 10, 10])
mol.set_calculator(calc)
# Check for correct input
if calc.dryrun_ok():
print('%s : %s ' % (mol.calc._label, mol.get_potential_energy()))
else:
print("Found error in input")
print(calc._error)
# Read all settings from previous calculation
mol = ase.io.castep.read_seed('%s/CO_LDA' % directory)
# Use the OTF pseudo-potential we have just generated
mol.calc.set_pspot('OTF')
# Change some settings
mol.calc.param.xc_functional = 'PBE'
# don't forget to set an appropriate label
mol.calc._label = 'CO_PBE'
# Recalculate the potential energy
print('%s : %s ' % (mol.calc._label, mol.get_potential_energy()))