All modern functional electronic devices are made from semiconductors. A semiconductor is a material with a moderate band gap. But to turn this semiconductor into a device, you need to be able to selectively inject and extract electrons from above and below this band gap, forming selective contacts.
The metrics we want to calculate are the Electron Affinity (energy gain by
adding an electron from infinite distance) and the Ionisation Potential (energy
required to remove an electron from the bulk to infinite distance).
In a metal there is no distinction between these values (the same energy state
is involved), and both are equivalent to the Work Function.
When calculating electronic structure, in a periodic system (where the electronic states are broadened bands) we talk of the Conduction Band Minimum (CBM) and Valence Band Maximum (VBM); in a molecular system (where for calculation we usually have an isolated molecule in vacuum) we talk in terms of the Lowest Unoccupied Molecular Orbital (LUMO) and Highest Occupied Molecular Orbital (HOMO).
The difference (in energy) between CBM and VBM, or similar LUMO and HOMO, gives us the Fundamental Gap. In a periodic system, this fundamental gap is a good estimate for the optical Band Gap of the material. This is because the absorption of a photon generates very loosely (negligibly) bound electron and hole states (the exciton). In a molecular system, this difference between the Fundamental and Band Gap (the exciton binding energy), can easily be 1 eV. The excited state (exciton) involves wave functions that are highly localised. The electron and hole are energetically bound; as there is large wavefunction overlap, the exchange interaction means that Singlet and Triplet states exist, and have different energies. Methods that explicitly treat the electron-hole interaction (such as Time-Dependent Density-Functional-Theory) are a much better estimate for the band gap in such systems.
But what are these orbital energies relative to? The naive expectation would be to the vacuum energy.
For a periodic calculation this can’t be, as there is no absolute reference in the periodic summations used to get the electrostatic potential. Often codes practically pin these values against the core states of the composite atoms used in the pseudo potential - so they ‘look’ about right, but woe betide anyone using these as actual work functions! However, you can define a ‘vacuum’ region in your simulation where no electrons are present and use the electrostatic potential in this region as the reference for your bulk values. My colleague Keith Butler has developed this method into a useful (and fun!) code named MacroDensity ([1] - https://github.com/WMD-group/MacroDensity).
With non-periodic (molecular) calculations, you do have an absolute reference - the vacuum is everywhere around you! But are the electron eigenvalues (i.e. HOMO / LUMO) actually what you want?
A key thing to remember is that the Density Functional Approximation is guaranteed to give you one thing - total energy. The Kohn-Sham eigenvectors [2] are auxiliary functions to enable the solving for the total energy of the system, in particular by providing orbitals which can be used to evaluate the kinetic energy of the electrons. The debate about whether the Kohn-Sham eigenvalues have any greater meaning has rumbled on for decades. It has become an almost ecumenical matter. It is usually referred to as “Koopmans’ Theorem” in the literature, which refers to the exact agreement in Hartree-Fock theory between the occupied orbitals and the ionisation energies (energetic cost to remove an electron). Depending on who you ask, this theory is correct for Hybrid Density Functional Theory for: all occupied and unoccupied orbitals, just the occupied orbitals, just the highest occupied orbital (HOMO).
Certainly in my testing with organic electron materials and hybrid density functionals, the HOMO often seems a good estimate for the ionisation potential. The LUMO is often not even correlated with variations in electron affinity (i.e. it is qualitatively and quantitatively wrong.)
So if we don’t trust our single-electron eigenvalues, what can we do?
A conceptually (and practically!) simple method is that of the ‘Delta SCF’.
The name is rather curious, but it seems to have stuck.
You calculate the total Self-Consistent-Field energy of a N electron system,
and deduct the total energy of a N+1 (or N-1) system.
This difference is the total energy upon adding (or removing) an electron, and
is therefore the Electron Affinity (or Ionisation Energy).
You do this at the geometry of the neutral molecule (apart from anything else,
this stops the difference in total energy requiring consideration of the
nuclear-nuclear repulsion which would change between two different geometries).
As ever in electronic structure theory this method ‘weighs the admiral by weighing the battleship with and without him on board’. It is thus extremely important that errors in the two total energies are carefully controlled.
An alternative, which chemists seem to particularly like, is to estimate the Ionisation Potential from the HOMO of a DFT calculation, and then add on the Band Gap of the material via a TDDFT calculation to get the Electron Affinity. It is somewhat surprising that this often gives a good agreement to the solid state energies. One possibility is that it represents some kind of ‘Pauling point’ (correct for the wrong reasons), that the exciton binding energy (and assorted systematic errors) in the TDDFT is somehow compensating for the delocalisation of the wavefunction in a solid, and the effect of the dielectric environment. This method I used in a 2009 study of various fullerene adduct isomers [3]. Later work (alas unpublished) on a much wider range of fullerenes suggested that the major variation in the Electron Affinity (~LUMO) was due to there being a constant Fundamental Gap between it and the HOMO.
Mark Casida (yes, that Casida!) has some interesting slides on why there should be such an agreement between TDDFT and Delta-SCF [4].
Doing it by hand with Gaussian
So you want to do a Delta-SCF calculation for a molecule-in-vacuum with Gaussian. You start with a carefully converged single point calculation:
#p sp b3lyp/6-31g* SCF(Tight,Conver=8) Integral(Grid=UltraFine)
Neutral calculation - for Delta-SCF.
0 1
A standard, single point, energy calculation, with upped convergence & the
two-electron integrals calculated on an ultrafine grid. The 0 1 is the
Charge Multiplicity. This is a Neutral (Charge=0) calculation in a spin Singlet
(Multiplicity=1) configuration.
We have made the convergence tighter (i.e. it will do more SCF cycles to get a lower noise answer) on an increased fidelity integration grid. See: http://www.gaussian.com/g_tech/g_ur/k_integral.htm
And now the same structure (no further nuclear relaxation), as a charged doublet.
If you want the energy of adding an electron to the ‘LUMO’ (the electron
affinity), you compare the neutral energy to the system with an added electron
(Charge=-1, the Anion). This is comparing N to N+1 electrons.
For a Delta-SCF calculation of the Ionisation Potential (IP), you compare to
the N-1 electron case, the cation.
#p sp b3lyp/6-31g* SCF(Tight,Conver=8) Integral(Grid=UltraFine)
Anion calculation - for Delta-SCF of the LUMO / Electron Affinity.
-1 2
And that is it. Now run these jobs, calculate the difference of the total SCF energy between them, and you have your estimate of the Ionisation Potential or the Electron Affinity.
#!/bin/sh
neutral=` grep "SCF Done" "${1}" | awk '{print $5}' `
charged=` grep "SCF Done" "${2}" | awk '{print $5}' `
echo -n "${1} "
echo -n "($neutral - $charged ) * 27.211" | bc -l
echo " eV"
Scripts for Gaussian
After doing this for the n’th time as part of this post, I thought I should sit down and write an all-singing script for automating the process.
The files are available as a Github Gist (https://gist.github.com/jarvist/616c8c7c3ea576d03d73).
From Gaussian .log outputs with optimised structures, generate .xyz files of
the coordinates with >extract_list_geoms.sh *.log.
Then generate your Delta-SCF .com Gaussian input files,
>g09_xyz_to_delta_SCF_coms_standalone.sh *.xyz.
Run the resulting .com files through Gaussian.
Analyse the outputs with > delta_SCF_extract_values.sh.
References
[1] Macrodensity - codes to read VASP LOCPOT files and calculate electrostatic potential derived quantities (https://github.com/WMD-group/MacroDensity)
[2] Kohn, Walter; Sham, Lu Jeu (1965). “Self-Consistent Equations Including Exchange and Correlation Effects”. Physical Review 140 (4A): A1133–A1138 (http://dx.doi.org/10.1103/PhysRev.140.A1133)
[3] Frost, J. M., Faist, M. A. and Nelson, J. (2010), Energetic Disorder in Higher Fullerene Adducts: A Quantum Chemical and Voltammetric Study. Adv. Mater., 22: 4881–4884. http://dx.doi.org/10.1002/adma.201002189
[4] Mark Casida - Reconciling Delta-SCF and TDDFT, Transparencies from a talk, 1999 (http://online.kitp.ucsb.edu/online/tddft-c99/casida/oh/01.html)
