Debugging

As a test, perform a three-dimensional calculation of a hydrogen atom as follows.

First, compute the potential energy of the periodic array of nuclei, proceeding as you did in Problem 1.1 with parameters sigma=0.25, S=[64;64;64], R=diag([16; 16; 16])), but placing just a single nucleus at the center of the cell. Be sure to subtract the self-energy, $1/(2\sqrt{\pi}\sigma)$.

Then, make a copy of ``dft.m'' called ``Hatom.m'' and replace the evaluation of Vdual there with the above octave expressions for rhohat (which represents J(rhonuc)) and Vdual for a Hydrogen nucleus (charge ZI=1) located at the center of the cell (rI=[8;8;8]) and run ``Hatom.m'' for a single state with one electron (gbl_f=1, Ns=1) with S=[64;64;64], R=diag([16; 16; 16]).

When your code functions properly (and you add the nuclear-nuclear potential energy), you should find a result within 1 millihartree of what you get either from your atomic code (with the VWN form for $\epsilon_{xc}$) or from looking directly on the NIST web site result. Now, though, you are doing a full three-dimensional calculation, without artificially imposing any spherical symmetry!