# pyscf: Why should the choice of pseudopotential determine the total energy, and not the physics?

I have limited experience with DFT, but as an exercise I have written a DFT program of an atom (from scratch using nothing more than python built in functions and numpy) by following Richard Martin's book on Electronic Structure calculations from 2004. Now I am trying to learn to use pyscf, but it has limited documentation.

• In the chapter of Richard Marin on pseudopotentials I read that --core states remain almost unchanged when a pseudopotential is introduced. So I do not anticipate that the total energy would change by a factor of 2 or 20 when starting with a pseudopotential.

However here is a script for bulk diamond (8 atoms per unit cell) after the atoms have been perturbed by a small amount and with lattice constant 3.57 using the Goedecker-Teter-Hutter (GTH) pseudopotential

import pyscf.pbc.dft as pbcdft
from pyscf.pbc import gto

cell = gto.Cell()
cell.a = '''
3.57  0       0
0       3.57  0
0       0       3.57'''
cell.atom =  '''C       -0.0059      -0.0384      -0.0080
C        0.8996       0.8834       0.8834
C        1.7947       1.7855       0.0279
C        2.6933       2.6673       0.8611
C        1.7998       0.0367       1.7496
C        2.7021       0.8854       2.6549
C       -0.0023       1.7795       1.8028
C        0.8828       2.6990       2.6548'''
cell.basis = 'gth-szv'
cell.pseudo = 'gth-lda'
cell.verbose = 4
cell.build()
mf=pbcdft.RKS(cell)
mf.kernel()


and using this script I calculate DFT energy (per unit cell) = -44.871672302677794. I assume that is -152.6275508125 eV per atom. But my reference says I should expect close to 9 eV at 5.687 Ångstrøm$$^3$$/atom. Can you see if it is it my reasoning that the pseudopotential shouldn't change the total energy by a factor of 17 that is wrong, or is it a problem with the PYSCF examples?

EDIT: I am assuming the input units of pyscf are in Ångstrøm, and the output is given in Hartree (because no units are quoted)

• If you want to calculate the energy per atom, you should divide the output energy by the number of atoms in your simulation (here, 8).
– user8097
Commented Jan 17 at 20:42
• Energies in DFT don't carry meaning on their own, but differences do, assuming what you're comparing was calculated with the same theory.
– user8097
Commented Jan 17 at 20:45
• Thanks. If there are 8 atoms, I divide by 8. But I also convert from hartree to eV. Commented Jan 17 at 21:26
• I cant to get the total energy. That's the trace of the product of the density matrix with the hamiltonian Commented Jan 17 at 21:26

The energy of a system, be it a molecule or a material, is dependent on the energy reference (or energy zero point) you choose. Common energy references are:

1. The state where all electrons (except for those implicitly described by pseudopotentials, if any) are ionized and have zero kinetic energy, and all nuclei are infinitely far from each other and also have zero kinetic energy. This convention is chosen by most molecular DFT codes.
2. The state where all atoms are free, neutral, and stationary. This convention is chosen by many periodic DFT codes. Some codes also assume the atoms are spin unpolarized and spherically symmetric.
3. The state where all atoms are in their reference elemental states (all carbon atoms as graphite, for example). This convention is chosen by many semiempirical codes.

Generally, the first convention yields very negative energies, because it takes a huge amount of energy to ionize all electrons of a system, or even just the valence electrons. A quick B3LYP/6-31G(d) calculation shows that the energy of a free carbon atom by this convention is about -38 Hartree when core electrons are treated explicitly, or about -5.5 Hartree when they are treated implicitly. The latter is quite close to your value (-153 eV), which proves that PySCF adopts the first convention even in periodic calculations.

By comparison, the second convention yields much less negative energies (but usually still negative, since most people study materials that are more stable than free atoms). The third convention generates values that are even smaller in magnitude, and they can often be positive.

Within a small enough community, it may be that all programs people use belong to one of the above categories, for example the second one. Then, within that community, it is unambiguous to talk about the absolute energy of a system (or at least its order of magnitude) without worrying about the energy reference. However, it is important to note that the other two conventions exist and are common in their respective communities. This is one of the reasons why energy differences are much more meaningful quantities to talk about in DFT calculations, compared to absolute energies.