My Community

When I am doing charged defect calculations and relaxing the ions, I would like to look at the total charge density (smoothed out a bit - singularities at the nuclei would not be good).
The ultimate goal is to use the total charge density in a Madelung-energy calculation to correct for the Coulomb interaction between the defect and it's unphysical lattice of copies (there are two recent papers on this by Freysoldt et al., but they suggest only using the electronic defect state charge).

I have adopted the following procedure, and I am wondering whether there might be some pitfalls with it, or whether there is a another way to do this. (Regarding the latter, there are lines one can uncomment in main.F that say they allow the total charge density to be written to the CHGCAR file, by adding the quantity DENCOR to the usual CHTOT. But my preliminary tests indicate this is not giving the true total charge; at least it differs from the results of the procedure I describe below, and it manifestly has the wrong total charge (it still gives the number of valence electrons)).

Here's what I do now:
I read the LOCPOT (generated with LVHAR=.TRUE.). I then FFT it. I take the Laplacian of it by multiplying by the (shifted) magnitude squared of the G vectors, with appropriate scaling, while at the same time multiplying each plane wave by an optional factor that decreases with G^2, to obtain a low pass filter.
I take the inverse FFT, and multiply by some carefully calculated scaling factors to obtain the total charge distribution. If I want, I can add a uniform charge to give the resulting distribution a net charge (nice for visualization, probably not wanted for Madelung type correction calculations). I feel I have been fairly careful in implementing this conversion, and have looked in the VASP code to make sure I know how the files are ordered and where the origin is and so on.

In principle this should get me the correct total charge distribution. I have the following questions
(1) Does this procedure suffer from serious problems resulting from the coarseness of the LOCPOT in the vicinity of the ions?
(2) In the limited testing I have done, the resulting shape of the charge distribution is qualitatively different from the result obtained by the uncommenting of the lines in main.F. (It's not just the total charge that is the problem, there are differences in the shape near the ions (ok - I tested it only on a single F- ion). My question is, does this indicate a problem in my LOCPOTtoCHG code? That is, is uncommenting those three lines really supposed to give the true total charge density (albeit off by a constant that is the number of valence electrons)? What exactly is supposed to physically be in the CHGCAR when those lines are uncommented?
(3) I use LOCPOT from LVHAR=T. Let me know if this does *not* include the (smoothed) electrostatic potential from the core electrons and the nucleus as well as the electron charge density. My looking at the manual and the code lead me to believe that it does include these things through the complicated internal PAW implementation, but I could very well be wrong.

Thank you very much,
David

From your post, it seems like you intend to find the charge density due to the valence electrons, core electrons, and the nucleus.

For calculations using PAW, the core charge density can be written out using the LAECHG tag. Note that for LVHAR, only the local part of the electrostatic potential is written out. That means the non-local part of the ionic potential will not be written in LOCPOT.

tlchan: Thank you for letting me know about LAECHG. However, I do not see how this can help me. With LAECHG there are three output files, two of which have a total charge close to the number of valence electrons while the third has a total charge which is, in my test on GaAs, larger than the total number of electrons (or protons) in the cell. So neither of these is what I want. I believe LAECHG is essentially providing a more accurate version the valence charge.

As tlchan correctly pointed out I do intend to find the *total* charge density, which is the charge of the valence electrons, the core electrons and the nucleus. It needs to be on a grid, and thus the nuclear and core electron charge must be spread out.

I believe the correct way to get the total charg is to take the Laplacian of the LOCPOT with LVHAR=TRUE. (LVTOT wouldn't make sense: I don't want exchange correlation potential, I want the electrostatic potential so I can use the microscopic Maxwell equation for rho and phi).

My primary question is simply, does LOCPOT actually contain a smoothed version of the true electrostatic potential?
The manual indicates this is so, but then the manual leaves out a detailed description of what is in CHGCAR and WAVECAR (that is, it doesn't explain that those aren't the true charges or wavefunctions). Do any VASP experts know the answer to my question about LOCPOT?

The answer to this question could help anyone seeking to correct for artificial electrostatic interaction between supercells for charged defect corrections. The suggested procedure by Freysoldt et al. (Phys. Status Solidi B, 248, 1067-1076 (2010)) uses the partial charge density of the "defect state". This only makes sense if the defect is very localized and clearly in the band gap, and if one does not relax the ions, which is often essential in a defect calculation. If one had access to the total charge, a Madelung-Freysoldt style correction could be implemented, I believe.

Thank you.

For LAECHG, you can do a search in the forum. I think people found that you need a very dense grid in order to reproduce the core charge accurately.

For the ionic potential, I doubt that LOCPOT will be helpful. In VASP, the ionic potential is represented by a pseudopotential, which is splitted into a local part and a non-local part. As noted in my post above, only the local part is plotted out in LOCPOT, the non-local part is an operator and cannot be plotted (the non-local part vanishes outside some cutoff radius. If you are only interested in the charge distribution outside the ionic core, then the LOCPOT is sufficient). So, I doubt that you can get the correct total charge by using the Poisson equation this way.
<span class='smallblacktext'>[ Edited Sat Apr 14 2012, 07:12AM ]</span>