How to obtain electron charge density per unit volume about a water molecule using SlowQuant?

Question

@Nike Dattani's answer to How to “get my feet wet” in Density Functional Theory by simulating a water molecule using Python recommends SlowQuant. I've used the install instructions for Linux on my MacOS and it seems to work.

Initially the tests failed because factorial2 has been moved from scipy.misc to scipy.special. I edited MIpython.py and BasisSet.py adding import scipy.special as spe and changing scp.factorial2 to spe.factorial2 in three places in each file and all tests now pass, and I get the following output:

$ python SlowQuant.py H2O.csv settingExample.csv

0.1089320182800293 INTEGRALS
Iter     Eel
1    -66.7536308450      -58.7512637840       6.71128773e+01      6.57412977e+00
2    -85.6091201524      -77.6067530914      -1.88554893e+01      6.02202187e+00      7.40515875e-01
3    -86.6677593084      -78.6653922474      -1.05863916e+00      1.07013952e-01      7.40515875e-01
0.006840229034423828 HF
4.0531158447265625e-06 PROPERTIES
2.1457672119140625e-06 Perturbation
2.1457672119140625e-06 QFIT
9.5367431640625e-07 CI
2.1457672119140625e-06 CC
0.1327359676361084 Execution time

So while I've been splashed my feet are not fully wet yet.

Question: What are the next steps to obtaining a 3D distribution of electronic charge in the molecule? This is all new to me.

I have seen Mulliken charges and Löwdin charges mentioned in the properties page of the documentation, perhaps something like this is what I'm after?

But when I added a line at the bottom of the script, I just got an error.

10;;;
1;1.63803684;1.136548823;0.0
8;0.0;-0.143225817;0.0
1;-1.63803684;1.136548823;0.0
Charge;Mulliken;;

I get

Traceback (most recent call last):
  File "SlowQuant.py", line 95, in <module>
    results = run(mol, set)
  File "SlowQuant.py", line 58, in run
    basis = BS.bassiset(input, set['basisset'])
  File "/your-name-here/SlowQuant/slowquant/basissets/BasisSet.py", line 195, in bassiset
    if basisload[j+1,0] == 'FOR' and writecheck == 1:
IndexError: index 176 is out of bounds for axis 0 with size 176

Ultimately I hope to draw some imaginary electric field lines myself and that's a challenge and topic of its own, but for right now I just need a way to get charge density per unit volume.

Answer 1

Firstly regarding the SlowQuant problems:

In the version of SlowQuant you are using, the program takes two files as input.

1. The molecules file with the coordinates, number of electrons, and atom IDs.
2. The settings file, see; https://slowquant.readthedocs.io/en/latest/Examples.html. It is in the settings file you would want to add "Charge;Mulliken"

If you have any problems with the code I can be reached by PM on https://www.reddit.com/user/ErikRK

Also if you want to add any features or have any bug fixed feel free to make a pull request on https://github.com/erikkjellgren/SlowQuant

Secondly you the charge density

The Mulliken charges and Lowdin charges is not quite what you are looking for.

If you want the charge density of the electrons in a point in space, this is what we also know as the electron density, and it is defined as follows:

$$\rho\left(\boldsymbol{r}\right)=\sum_{pq}\Omega_{pq}\left(\boldsymbol{r}\right)D_{pq}$$

Here $D$ is the density matrix and $\Omega$ is the overlap distribution of two atomic orbitals and given as:

$$\Omega_{pq}(\boldsymbol{r}) = \phi_p(\boldsymbol{r})\phi_q(\boldsymbol{r})$$

You will therefore have to evaluate the value of the basis functions in a point in space, and then weight them by the molecular orbital coefficients.