Usage in ASE#

The Python API of tblite natively support integration with the atomic simulation environment (ASE). By constructing a calculator most functionality of ASE is readily available. For details on building the Python API checkout the installation guide.

ASE calculator#

class tblite.ase.TBLite(*args: Any, **kwargs: Any)[source]#

ASE calculator for using xTB Hamiltonians from the tblite library. Supported properties by this calculator are:

  • energy (free_energy)

  • forces

  • stress

  • dipole

  • charges

Supported keywords are

Keyword

Default

Description

method

“GFN2-xTB”

Underlying method for energy and forces

charge

None

Total charge of the system

multiplicity

None

Total multiplicity of the system

accuracy

1.0

Numerical accuracy of the calculation

electronic_temperature

300.0

Electronic temperatur in Kelvin

max_iterations

250

Iterations for self-consistent evaluation

initial_guess

“sad”

Initial guess for wavefunction (sad or eeq)

mixer_damping

0.4

Damping parameter for self-consistent mixer

electric_field

None

Uniform electric field vector (in V/A)

spin_polarization

None

Spin polarization (scaling factor)

cache_api

True

Reuse generate API objects (recommended)

verbosity

1

Set verbosity of printout

Example

An ASE calculator can be constructed by using the TBLite class provided by the tblite.ase module. For example to perform a single point calculation for a CO2 crystal use

>>> from tblite.ase import TBLite
>>> from ase.atoms import Atoms
>>> import numpy as np
>>> atoms = Atoms(
...     symbols="C4O8",
...     positions=np.array(
...         [
...             [0.9441259872, 0.9437851680, 0.9543505632],
...             [3.7179966528, 0.9556570368, 3.7316862240],
...             [3.7159517376, 3.7149292800, 0.9692330016],
...             [0.9529872864, 3.7220864832, 3.7296981120],
...             [1.6213905408, 1.6190616096, 1.6313879040],
...             [0.2656685664, 0.2694175776, 0.2776540416],
...             [4.3914553920, 1.6346256864, 3.0545920000],
...             [3.0440834880, 0.2764611744, 4.4080419264],
...             [4.3910577696, 3.0416409504, 0.2881058304],
...             [3.0399936576, 4.3879335936, 1.6497353376],
...             [0.2741322432, 4.4003734944, 3.0573754368],
...             [1.6312174944, 3.0434586528, 4.4023048032],
...         ]
...     ),
...     cell=np.array([5.68032, 5.68032, 5.68032]),
...     pbc=np.array([True, True, True]),
... )
>>> atoms.calc = TBLite(method="GFN1-xTB")
>>> atoms.get_potential_energy()  # result in eV
-1257.0943962462964

The resulting calculator can be used like most ASE calculator, e.g. for optimizing geometries.

calculate(atoms: ase.atoms.Atoms | None = None, properties: List[str] = None, system_changes: List[str] = ase.calculators.calculator.all_changes) None[source]#

Perform actual calculation with by calling the TBLite API

Example

>>> from ase.build import molecule
>>> from tblite.ase import TBLite
>>> calc = TBLite(method="GFN2-xTB")
>>> calc.calculate(molecule("H2O"))
>>> calc.get_potential_energy()
-137.96777625229421
>>> calc.calculate(molecule("CH4"))
>>> calc.get_potential_energy()
-113.60956621093894
Raises:
  • ase.calculators.calculator.InputError – on invalid input passed to the interface module

  • ase.calculators.calculator.CalculationFailed – in case of an RuntimeError in the library

default_parameters: Dict[str, Any] = {'accuracy': 1.0, 'cache_api': True, 'charge': None, 'electric_field': None, 'electronic_temperature': 300.0, 'guess': 'sad', 'max_iterations': 250, 'method': 'GFN2-xTB', 'mixer_damping': 0.4, 'multiplicity': None, 'spin_polarization': None, 'verbosity': 1}#

Default parameters

implemented_properties: List[str] = ['energy', 'forces', 'charges', 'dipole', 'stress']#

Properties calculator can handle (energy, forces, …)

reset() None[source]#

Clear all information from old calculation. This will only remove the cached API objects in case the cache_api is set to False.

set(**kwargs) dict[source]#

Set new parameters to TBLite. Will automatically reconstruct the underlying model in case critical parameters change.

Example

>>> from ase.build import molecule
>>> from tblite.ase import TBLite
>>> atoms = molecule("H2O")
>>> atoms.calc = TBLite(method="GFN2-xTB")
>>> atoms.get_potential_energy()
-137.96777625229421
>>> atoms.calc.set(method="GFN1-xTB")
{'method': 'GFN1-xTB'}
>>> atoms.get_potential_energy()
-156.9675057724589