ASE Support#
ASE calculator implementation
for the dftd4
program.
This module provides a basic single point calculator implementations
to integrate the dftd4
API into existing ASE workflows.
To use DFTD4 as dispersion correction the ase.calculators.mixing
module can be used to combine DFTD4 with a DFT calculator using
the SumCalculator
.
Supported properties by this calculator are:
energy (free_energy)
forces
stress
Supported keywords are
Keyword |
Default |
Description |
---|---|---|
method |
None |
Method to calculate dispersion for |
params_tweaks |
None |
Optional dict with the damping parameters |
cache_api |
True |
Reuse generate API objects (recommended) |
model |
d4 |
Used dispersion Model (D4S or D4 (default)) |
The params_tweaks dict contains the damping parameters, at least s8, a1 and a2 must be provided
Tweakable parameter |
Default |
Description |
---|---|---|
s6 |
1.0 |
Scaling of the dipole-dipole dispersion |
s8 |
None |
Scaling of the dipole-quadrupole dispersion |
s9 |
1.0 |
Scaling of the three-body dispersion energy |
a1 |
None |
Scaling of the critical radii |
a2 |
None |
Offset of the critical radii |
alp |
16.0 |
Exponent of the zero damping (ATM only) |
Either method or s8, a1 and a2 must be provided, s9 can be used to overwrite the ATM scaling if the method is provided in the model. Disabling the three-body dispersion (s9=0.0) changes the internal selection rules for damping parameters of a given method and prefers special two-body only damping parameters if available!
Example
>>> from ase.build import molecule
>>> from dftd4.ase import DFTD4
>>> atoms = molecule('H2O')
>>> atoms.calc = DFTD4(method="TPSS")
>>> atoms.get_potential_energy()
-0.007310393443152083
>>> atoms.calc.set(method="PBE")
{'method': 'PBE'}
>>> atoms.get_potential_energy()
-0.005358475432239303
>>> atoms.get_forces()
array([[-0. , -0. , 0.00296845],
[-0. , 0.00119152, -0.00148423],
[-0. , -0.00119152, -0.00148423]])
- class dftd4.ase.DFTD4(*args, **kwargs)[source]#
ASE calculator for DFT-D4 related methods. The DFTD4 class can access all methods exposed by the
dftd4
API.Example
>>> from ase.build import molecule >>> from ase.calculators.mixing import SumCalculator >>> from ase.calculators.nwchem import NWChem >>> from dftd4.ase import DFTD4 >>> atoms = molecule('H2O') >>> atoms.calc = SumCalculator([DFTD4(method="PBE"), NWChem(xc="PBE")])
- add_calculator(other)[source]#
Convenience function to allow DFTD4 to combine itself with another calculator by returning a SumCalculator:
Example
>>> from ase.build import molecule >>> from ase.calculators.emt import EMT >>> from dftd4.ase import DFTD4 >>> atoms = molecule("C60") >>> atoms.calc = DFTD4(method="pbe").add_calculator(EMT()) >>> atoms.get_potential_energy() 6.348142387048062 >>> [calc.get_potential_energy() for calc in atoms.calc.calcs] [-6.015477436263984, 12.363619823312046]