Hello AJAK,
Your input file is correct, but there are a number of parameters and some features of the DFT-LDA approximations that explain why the answer is not the exact value (-0.5). This calculation uses the default option for the "xc" variable, i.e. xc=LDA. In the local density approximation (LDA) the exchange and correlation energy is described approximately using an expression that depends only on the electronic charge density locally. This means that the one-electron effective potential felt by an electron depends on the electronic charge density of all electrons. In the hydrogen atom (with only one electron), this leads to a spurious "self-interaction", i.e. the effective potential seen by the only electron depends on its own charge, i.e. it is interacting with itself. This is obviously a bad approximation for the hydrogen atom, but turns out to be a more reasonable approximation in larger, many-electron systems.
A number of modifications to the input file can be made to improve the result, and I describe them below:
1) The value of the ecut variable (plane wave energy cutoff) is quite low. It is important to check that the results are converged with respect to the plane wave energy cutoff. For hydrogen, a value of 35 Rydberg is usually adequate:
Code: Select all
set cell 20 0 0 0 20 0 0 0 20
species hydrogen hydrogen.xml
atom H1 hydrogen 10 10 10
set ecut 35
set wf_dyn PSDA
set ecutprec 4
randomize_wf
set scf_tol 1.e-8
run 0 200
(note also the use of the "scf_tol" variable which will stop the calculation as soon as a tolerance of 1.e-8 is reached).
This calculation yields etotal = -0.44297092.
2) The default option in Qbox is to ignore the spin of electrons, and describe orbitals without a spin degree of freedom. This is adequate in many situations (such as closed-shell systems), but is definitely bad for hydrogen (only one electron). This can be changed by setting the variable nspin to 2
Code: Select all
set cell 20 0 0 0 20 0 0 0 20
species hydrogen hydrogen.xml
atom H1 hydrogen 10 10 10
set ecut 35
set nspin 2
set wf_dyn PSDA
set ecutprec 4
randomize_wf
set scf_tol 1.e-8
run 0 200
This yields etotal=-0.47555553, which is better, but the self-interaction error of LDA is still sizeable.
3) If we replace the LDA by the Hartree-Fock approximation, it is much better for hydrogen since it is free of self-interaction errors. Indeed, for hydrogen it should give the exact result:
Code: Select all
set cell 20 0 0 0 20 0 0 0 20
species hydrogen hydrogen.xml
atom H1 hydrogen 10 10 10
set ecut 35
set nspin 2
set xc HF
set wf_dyn PSDA
set ecutprec 4
randomize_wf
set scf_tol 1.e-8
run 0 200
This calculation yields etotal=-0.49744650. This still differs from the exact answer, but now this is due to remaining numerical approximations. First, increasing the value of ecut to 50 Ry yields etotal=-0.49954032. Increasing it further to 70 Ry yields etotal=-0.50043248.
We now have a value of etotal which is
below the exact value. This is due to the fact that we are dealing with a periodic unit cell. We have assumed that the hydrogen atom is well isolated from periodic replicas by setting the cell size to 20 Bohr. There is still a very weak binding interaction between hydrogen atoms even at that distance. As a result, etotal is slightly below -0.5. We can test that further by increasing the cell size to 30 Bohr:
Code: Select all
set cell 30 0 0 0 30 0 0 0 30
species hydrogen hydrogen.xml
atom H1 hydrogen 10 10 10
set ecut 70
set nspin 2
set xc HF
set wf_dyn PSDA
set ecutprec 4
randomize_wf
set scf_tol 1.e-8
run 0 200
This calculation yields etotal=-0.49987738.
In summary, it is important to choose carefully the level of approximation of the exchange-correlation energy (LDA, or HF, or other), to check the accuracy of the plane wave basis by increasing the value of ecut, and to check the effect of periodic boundary conditions by increasing the unit cell size.
Hope this helps,
Francois