Calculation of atomization energies¶
Warning: mainstream DFT is unable to describe correctly the electronic states of isolated atoms (especially of transition metals). See http://dx.doi.org/10.1063/1.2723118 . This usually manifests itself as SCF convergence problems. Please consult literature before reporting such problems on the mailing lists.
The following script will calculate the atomization energy of a hydrogen molecule:
# Creates: atomization.txt from __future__ import print_function from ase import Atoms from ase.parallel import paropen as open from gpaw import GPAW, PW, FermiDirac a = 10. # Size of unit cell (Angstrom) c = a / 2 # Hydrogen atom: atom = Atoms('H', positions=[(c, c, c)], magmoms=, cell=(a, a + 0.0001, a + 0.0002)) # Break cell symmetry # gpaw calculator: calc = GPAW(mode=PW(), xc='PBE', hund=True, eigensolver='rmm-diis', # This solver can parallelize over bands occupations=FermiDirac(0.0, fixmagmom=True), txt='H.out', ) atom.set_calculator(calc) e1 = atom.get_potential_energy() calc.write('H.gpw') # Hydrogen molecule: d = 0.74 # Experimental bond length molecule = Atoms('H2', positions=([c - d / 2, c, c], [c + d / 2, c, c]), cell=(a, a, a)) calc.set(txt='H2.out') calc.set(hund=False) # No hund rule for molecules molecule.set_calculator(calc) e2 = molecule.get_potential_energy() calc.write('H2.gpw') fd = open('atomization.txt', 'w') print(' hydrogen atom energy: %5.2f eV' % e1, file=fd) print(' hydrogen molecule energy: %5.2f eV' % e2, file=fd) print(' atomization energy: %5.2f eV' % (2 * e1 - e2), file=fd) fd.close()
Atoms object containing one hydrogen atom with a
magnetic moment of one, is created. Next, a GPAW calculator is
created. The calculator will do a calculation
using the PBE exchange-correlation functional, and output from the
calculation will be written to a file
H.out. The calculator is
hooked on to the
atom object, and the energy is calculated (the
e1 variable). Finally, the result of the calculation
(wavefunctions, densities, …) is saved to a file.
molecule object is defined, holding the hydrogen molecule at
the experimental lattice constant. The calculator used for the atom
calculation is used again for the molecule caluclation - only the
filename for the output file needs to be changed to
extract the energy into the
If we run this script, we get the following result:
hydrogen atom energy: -1.07 eV hydrogen molecule energy: -6.64 eV atomization energy: 4.50 eV
According to Blaha et al. 1, an all-electron calculation with PBE gives an atomization energy of 4.54 eV, which is in perfect agreement with our PAW result.
The energy of the spin polarized hydrogen atom is -1.09 eV. If we do
the calculation for the atom with
hund=False, then we get
almost 0 eV. This number should converge to exactly zero for a very
large cell and a very high grid-point density, because the energy of a
non spin-polarized hydrogen atom is the reference energy.
S. Kurth, J. P. Perdew, and P. Blaha, Int. J. Quantum Chem. 75, 889 (1999)