# -*- coding: utf-8 -*-
"""Contains a class for DFTB calculations using DFTBaby."""
import numpy as np
from .iface import iface
from .. import atoms
[docs]class iface_dftbaby(iface):
"""
Interface to DFTBaby for DFTB calculations, inherits from `metafalcon.interfaces.iface.iface`.
Parameters
----------
symbols : list of str
atomic symbols of the molecule
coords : np.2darray (shape: (N, 3))
atomic xyz-coordinates of the molecule (bohr)
**kwargs :
See below
Keyword Arguments
-----------------
charge : int
atomic charge
Attributes
----------
symbols : list of str
atomic symbols of the molecule
coords : np.2darray (shape: (N, 3))
atomic xyz-coordinates of the molecule (bohr)
energy : float
electronic energy of the molecule calculated most recently (hartrees)
gradient : np.2darray (shape: (N, 3))
electronic energy gradient of the molecule with current coordinates (hartree/bohr)
charge : int
atomic charge
dftbpes : DFTB.PES.PotentialEnergySurfaces
DFTBaby PotentialEnergySurfaces object for the calculation of energies and gradients
See Also
--------
iface.set_symbols()
iface.set_coords()
iface.get_energy()
iface.get_gradient()
"""
def __init__(self, symbols, coords, **kwargs):
"""Construct DFTBaby interface object."""
from DFTB.PES import PotentialEnergySurfaces
iface.__init__(self, symbols, coords)
self.charge = kwargs["charge"]
# initialize dftbaby object
atomlist = [(atoms.elements[sym.capitalize()], (0., 0., 0.)) for sym in self.symbols]
self.dftbpes = PotentialEnergySurfaces(atomlist, charge=self.charge)
# dftbaby needs one excited states calculation to set all variables
self.dftbpes.getEnergies(self.coords.flatten())
[docs] def calculate(self):
"""Calculate energies and gradients for the current set of coordinates."""
from DFTB.DFTB2 import SCFNotConverged
try:
self.energy = self.dftbpes.getEnergy_S0(self.coords.flatten())
# save CPKS intermediates during gradient calculation for possible Mullken charge CVs
gradVrep, gradE0, gradExc = self.dftbpes.grads.gradient(I=0, save_intermediates_CPKS=1)
self.gradient = (gradVrep + gradE0 + gradExc).reshape(self.coords.shape)
self.converged = True
except SCFNotConverged as e:
self.energy = 0.0
self.gradient = np.zeros_like(self.coords)
self.converged = False
print(e)