API Documentation¶
Important
All properties exchanged with the xtb
API are given in atomic units.
For integrations with other frameworks the unit conventions might differ and require conversion.
Contents
Calculation Environment¶
-
class
xtb.interface.
Environment
[source]¶ Wraps an API object representing a TEnvironment class in
xtb
. The API object is constructed automatically and deconstructed on garbage collection, it stores the IO configuration and the error log of the API.All API calls require an environment object, usually this is done automatically as all other classes inherent from the calculation environment.
Example
>>> from xtb.libxtb import VERBOSITY_FULL >>> from xtb.interface import Environment >>> env = Environment() >>> env.set_output("error.log") >>> env.set_verbosity(VERBOSITY_FULL) >>> if env.check != 0: ... env.show("Error message") ... >>> env.release_output()
-
check
() → int[source]¶ Check current status of calculation environment
Example
>>> if env.check() != 0: ... raise XTBException("Error occured in the API")
-
get_error
(message: Optional[str] = None) → str[source]¶ Check for error messages
Example
>>> if env.check() != 0: ... raise XTBException(env.get_error())
-
Molecular Structure Data¶
-
class
xtb.interface.
Molecule
(numbers, positions, charge: Optional[float] = None, uhf: Optional[int] = None, lattice=None, periodic=None)[source]¶ Represents a wrapped TMolecule API object in
xtb
. The molecular structure data object has a fixed number of atoms and immutable atomic identifiers.Example
>>> from xtb.interface import Molecule >>> import numpy as np >>> numbers = np.array([8, 1, 1]) >>> positions = np.array([ ... [ 0.00000000000000, 0.00000000000000,-0.73578586109551], ... [ 1.44183152868459, 0.00000000000000, 0.36789293054775], ... [-1.44183152868459, 0.00000000000000, 0.36789293054775]]) ... >>> mol = Molecule(numbers, positions) >>> len(mol) 3 >>> mol.update(np.zeros((len(mol), 3))) # will fail nuclear fusion check xtb.interface.XTBException: Update of molecular structure failed: -1- xtb_api_updateMolecule: Could not update molecular structure >>> mol.update(positions)
Raises: ValueError
– on invalid input on the Python side of the APIXTBException
– on errors returned from the API
-
update
(positions: numpy.ndarray, lattice: Optional[numpy.ndarray] = None) → None[source]¶ Update coordinates and lattice parameters, both provided in atomic units (Bohr). The lattice update is optional also for periodic structures.
Generally, only the cartesian coordinates and the lattice parameters can be updated, every other modification, regarding total charge, total spin, boundary condition, atomic types or number of atoms requires the complete reconstruction of the object.
Raises: ValueError
– on invalid input on the Python side of the APIXTBException
– on errors returned from the API, usually from nuclear fusion check
Single Point Calculator¶
-
class
xtb.interface.
Calculator
(param: xtb.interface.Param, numbers: List[int], positions: List[float], charge: Optional[float] = None, uhf: Optional[int] = None, lattice: Optional[List[float]] = None, periodic: Optional[List[bool]] = None)[source]¶ This calculator represents a calculator object in the
xtb
API and provides access to all methods implemented with a unified interface. The API object must be loaded with a parametrisation before it can be used in any other API request.The parametrisation loading is included in the initialization in this class, which has the advantage that all API functionality is readily available, the downside is that a calculator object on the Python side can only carry one distinct parametrisation, which is not allowed to change.
Examples
>>> from xtb.libxtb import VERBOSITY_MINIMAL >>> from xtb.interface import Calculator, Param >>> import numpy as np >>> numbers = np.array([8, 1, 1]) >>> positions = np.array([ ... [ 0.00000000000000, 0.00000000000000,-0.73578586109551], ... [ 1.44183152868459, 0.00000000000000, 0.36789293054775], ... [-1.44183152868459, 0.00000000000000, 0.36789293054775]]) ... >>> calc = Calculator(Param.GFN2xTB, numbers, positions) >>> calc.set_verbosity(VERBOSITY_MINIMAL) >>> res = calc.singlepoint() # energy printed is only the electronic part 1 -5.1027888 -0.510279E+01 0.421E+00 14.83 0.0 T 2 -5.1040645 -0.127572E-02 0.242E+00 14.55 1.0 T 3 -5.1042978 -0.233350E-03 0.381E-01 14.33 1.0 T 4 -5.1043581 -0.602769E-04 0.885E-02 14.48 1.0 T 5 -5.1043609 -0.280751E-05 0.566E-02 14.43 1.0 T 6 -5.1043628 -0.188160E-05 0.131E-03 14.45 44.1 T 7 -5.1043628 -0.455326E-09 0.978E-04 14.45 59.1 T 8 -5.1043628 -0.572169E-09 0.192E-05 14.45 3009.1 T SCC iter. ... 0 min, 0.022 sec gradient ... 0 min, 0.000 sec >>> res.get_energy() -5.070451354836705 >>> res.get_gradient() [[ 6.24500451e-17 -3.47909735e-17 -5.07156941e-03] [-1.24839222e-03 2.43536791e-17 2.53578470e-03] [ 1.24839222e-03 1.04372944e-17 2.53578470e-03]]
Raises: XTBException
– on errors encountered in API or while performing calculations-
set_accuracy
(accuracy: float) → None[source]¶ Set numerical accuracy for calculation, ranges from 1000 to 0.0001, values outside this range will be cutted with warning placed in the error log, which can be retrieved by get_error() but will not trigger check().
Example
>>> calc.set_accuracy(1.0)
-
set_electronic_temperature
(etemp: int) → None[source]¶ Set electronic temperature in K for tight binding Hamiltonians, values smaller or equal to zero will be silently ignored by the API.
Example
>>> calc.set_electronic_temperature(300.0)
-
set_external_charges
(numbers: numpy.ndarray, charges: numpy.ndarray, positions: numpy.ndarray) → None[source]¶ Set an external point charge field
-
set_max_iterations
(maxiter: int) → None[source]¶ Set maximum number of iterations for self-consistent charge methods, values smaller than one will be silently ignored by the API. Failing to converge in a given number of cycles is not necessarily reported as an error by the API.
Example
>>> calc.set_max_iterations(100)
-
set_solvent
(solvent: Optional[xtb.interface.Solvent] = None) → None[source]¶ Add/Remove a solvation model to/from calculator
Example
>>> from xtb.utils import get_solvent, Solvent ... >>> calc.set_solvent(Solvent.h2o) # Set solvent to water with enumerator >>> calc.set_solvent() # Release solvent again >>> calc.set_solvent(get_solvent("CHCl3")) # Find correct enumerator
-
singlepoint
(res: Optional[xtb.interface.Results] = None, copy: bool = False) → xtb.interface.Results[source]¶ Perform singlepoint calculation, note that the a previous result is overwritten by default.
Example
>>> res = calc.singlepoint() >>> res = calc.singlepoint(res) >>> calc.singlepoint(res) # equivalent to the above >>> new = calc.singlepoint(res, copy=True)
-
Calculation Results¶
-
class
xtb.interface.
Results
(res: Union[xtb.interface.Molecule, Results])[source]¶ Holds
xtb
API object containing results from a single point calculation. It can be queried for indiviual properties or used to restart calculations. Note that results from different methods are generally incompatible, the API tries to be as clever as possible about this and will usually automatically reallocate missmatched results objects as necessary.The results objects is connected to its own, independent environment, giving it its own error stack and IO infrastructure.
Example
>>> from xtb.libxtb import VERBOSITY_MINIMAL >>> from xtb.interface import Calculator, Param >>> import numpy as np >>> numbers = np.array([8, 1, 1]) >>> positions = np.array([ ... [ 0.00000000000000, 0.00000000000000,-0.73578586109551], ... [ 1.44183152868459, 0.00000000000000, 0.36789293054775], ... [-1.44183152868459, 0.00000000000000, 0.36789293054775]]) ... >>> calc = Calculator(Param.GFN2xTB, numbers, positions) >>> calc.set_verbosity(VERBOSITY_MINIMAL) >>> res = calc.singlepoint() # energy printed is only the electronic part 1 -5.1027888 -0.510279E+01 0.421E+00 14.83 0.0 T 2 -5.1040645 -0.127572E-02 0.242E+00 14.55 1.0 T 3 -5.1042978 -0.233350E-03 0.381E-01 14.33 1.0 T 4 -5.1043581 -0.602769E-04 0.885E-02 14.48 1.0 T 5 -5.1043609 -0.280751E-05 0.566E-02 14.43 1.0 T 6 -5.1043628 -0.188160E-05 0.131E-03 14.45 44.1 T 7 -5.1043628 -0.455326E-09 0.978E-04 14.45 59.1 T 8 -5.1043628 -0.572169E-09 0.192E-05 14.45 3009.1 T SCC iter. ... 0 min, 0.022 sec gradient ... 0 min, 0.000 sec >>> res.get_energy() -5.070451354836705 >>> res.get_gradient() [[ 6.24500451e-17 -3.47909735e-17 -5.07156941e-03] [-1.24839222e-03 2.43536791e-17 2.53578470e-03] [ 1.24839222e-03 1.04372944e-17 2.53578470e-03]] >>> res = calc.singlepoint(res) 1 -5.1043628 -0.510436E+01 0.898E-08 14.45 0.0 T 2 -5.1043628 -0.266454E-14 0.436E-08 14.45 100000.0 T 3 -5.1043628 0.177636E-14 0.137E-08 14.45 100000.0 T SCC iter. ... 0 min, 0.001 sec gradient ... 0 min, 0.000 sec >>> res.get_charges() [-0.56317912 0.28158956 0.28158956]
Raises: XTBException
– in case the requested property is not present in the results object-
get_bond_orders
() → numpy.ndarray[source]¶ Query singlepoint results object for bond orders
Example
>>> res.get_bond_orders() [[0.00000000e+00 9.20433501e-01 9.20433501e-01] [9.20433501e-01 0.00000000e+00 2.74039053e-04] [9.20433501e-01 2.74039053e-04 0.00000000e+00]]
-
get_charges
() → numpy.ndarray[source]¶ Query singlepoint results object for partial charges in e
Example
>>> get_charges() [-0.56317913 0.28158957 0.28158957]
-
get_dipole
() → numpy.ndarray[source]¶ Query singlepoint results object for dipole in e·Bohr
Example
>>> get_dipole() [-4.44089210e-16 1.44419023e-16 8.89047667e-01]
-
get_energy
() → float[source]¶ Query singlepoint results object for energy in Hartree
Example
>>> res.get_energy() -5.070451354836705
-
get_gradient
() → numpy.ndarray[source]¶ Query singlepoint results object for gradient in Hartree/Bohr
Example
>>> res.get_gradient() [[ 6.24500451e-17 -3.47909735e-17 -5.07156941e-03] [-1.24839222e-03 2.43536791e-17 2.53578470e-03] [ 1.24839222e-03 1.04372944e-17 2.53578470e-03]]
-
get_number_of_orbitals
() → int[source]¶ Query singlepoint results object for the number of basis functions
Example
>>> res.get_number_of_orbitals() 6
-
get_orbital_coefficients
() → numpy.ndarray[source]¶ Query singlepoint results object for orbital coefficients
Example
>>> res.get_orbital_coefficients() array([[-7.94626768e-01, 6.38378239e-16, 4.52990407e-01, -6.38746369e-16, -8.35495085e-01, -4.44089210e-16], [ 2.77555756e-17, -6.97332245e-01, 7.49400542e-16, 1.88136491e-17, 7.21644966e-16, -9.60006511e-01], [ 2.17336312e-16, -1.08051945e-16, -1.11598977e-15, -1.00000000e+00, 5.74153329e-17, 3.30330107e-17], [-8.67578876e-02, -9.71445147e-16, -8.05763104e-01, 7.71702239e-16, -7.18690020e-01, -4.71844785e-16], [-1.84540457e-01, -3.54572323e-01, -2.39090946e-01, 2.87533552e-16, 7.68757806e-01, 9.02845514e-01], [-1.84540457e-01, 3.54572323e-01, -2.39090946e-01, 2.01021058e-16, 7.68757806e-01, -9.02845514e-01]])
-
get_orbital_eigenvalues
() → numpy.ndarray[source]¶ Query singlepoint results object for orbital energies in Hartree
Example
>>> res.get_orbital_eigenvalues() array([-0.68087967, -0.56667693, -0.51373083, -0.44710101, 0.08394016, 0.24142397])
-
Available Calculation Methods¶
-
class
xtb.interface.
Param
[source]¶ Possible parametrisations for the Calculator class
-
GFN0xTB
= 3¶ Experimental non-self-consistent extended tight binding Hamiltonian using classical electronegativity equilibration electrostatics and extended Hückel Hamiltonian.
Geometry, frequency and non-covalent interactions parametrisation for elements up to Z=86.
Requires the param_gfn0-xtb.txt parameter file in the
XTBPATH
environment variable to load!See: P. Pracht, E. Caldeweyher, S. Ehlert, S. Grimme, ChemRxiv, 2019, preprint. DOI: 10.26434/chemrxiv.8326202.v1
-
GFN1xTB
= 2¶ Self-consistent extended tight binding Hamiltonian with isotropic second order electrostatic contributions and third order on-site contributions.
Geometry, frequency and non-covalent interactions parametrisation for elements up to Z=86.
Cite as: S. Grimme, C. Bannwarth, P. Shushkov, J. Chem. Theory Comput., 2017, 13, 1989-2009. DOI: 10.1021/acs.jctc.7b00118
-
GFN2xTB
= 1¶ Self-consistent extended tight binding Hamiltonian with anisotropic second order electrostatic contributions, third order on-site contributions and self-consistent D4 dispersion.
Geometry, frequency and non-covalent interactions parametrisation for elements up to Z=86.
Cite as: C. Bannwarth, S. Ehlert and S. Grimme., J. Chem. Theory Comput., 2019, 15, 1652-1671. DOI: 10.1021/acs.jctc.8b01176
-
GFNFF
= 4¶ General force field parametrized for geometry, frequency and non-covalent interactions up to Z=86.
xtb
API support is currently experimental.Cite as: S. Spicher and S. Grimme, Angew. Chem. Int. Ed., 2020, 59, 15665–15673. DOI: 10.1002/anie.202004239
-
IPEAxTB
= 5¶ Special parametrisation for the GFN1-xTB Hamiltonian to improve the description of vertical ionisation potentials and electron affinities. Uses additional diffuse s-functions on light main group elements. Parametrised up to Z=86.
Cite as: V. Asgeirsson, C. Bauer and S. Grimme, Chem. Sci., 2017, 8, 4879. DOI: `10.1039/c7sc00601b <https://dx.doi.org/10.1039/c7sc00601b`_
-
-
utils.
get_method
() → Optional[xtb.interface.Param]¶ Return the correct parameter enumerator for a string input.
Example
>>> get_method('GFN2-xTB') <Param.GFN2xTB: 1> >>> get_method('gfn2xtb') <Param.GFN2xTB: 1> >>> get_method('GFN-xTB') is None True >>> get_method('GFN1-xTB') is None False