Add Hirshfeld charge analysis #494
Conversation
charge (useful for debugging)
There was a problem hiding this comment.
This looks great! I think we can put this radial interpolator to more use, e.g. for Susi's Superposition of Atomic Potentials (SAP) initial guess. And these LDA atomic densities are probably better than the 3-21G densities used in the current SAD initial guess, so I guess we can replace those 🙂
|
Are the SAD input densities also spherically symmetric? If you think the radial interpolator is useful for other people too, we should probably put it in a different directory do you have some suggestions? |
|
Currently the SAD densities are not symmetric, because they are precomputed using unrestricted LDA in a tiny 3-21G basis, but I think it would be better if they were symmetric. Perhaps move it to a folder under |
|
Yes, in that case I think it is better when they are symmetric. I can polish the interpolator a bit, implement some more interpolation modes and move it to utils. |
|
I moved specific parts of my implementation into more accessible locations:
I also added another extrapolation mode in the poynomial interpolator that is useful for interpolating densities directly: No extrapolation to the left and returning 0 if x > xmax. The building blocks to replace the SAD densities by our new densities should now be in place and it should be quite easy to implement it in the future |
I implemented Hirshfeld charge partitioning. A method to compute partial charges of atoms in molecules which is implemented in many electronic structure code. The Hirshfeld charge of atom$i$ is defined as:
$$q_i=Z_i-\int\frac{\rho^0_i(\mathbf{r})}{\sum_j\rho^0_j(\mathbf{r})}\rho(r)d\mathbf{r}$$ $Z_i$ is the atomic number of element $i$ , $\rho$ is the molecular density and $\rho_i^0$ is density of $i$ as an isolated atom.
where
The spherically symmetric reference charge densities were obtained with an LDA calculation and are tabulated in share/hirshfeld/lda/{element_name}.density
To improve the numerical stability, the logarithm of the density is interpolated internally and the logsumexp trick was used to get rid of 0/0 issues in regions far away from all nuclei in the function in the integral above.