shankar1729 / jdftx

JDFTx: software for joint density functional theory
http://jdftx.org
81 stars 54 forks source link

Coulomb cutoff implementation in JDFTx #188

Closed MSJavaScript closed 2 years ago

MSJavaScript commented 2 years ago

How coulomb cutoff in JDFTx works ? Can you give me some reference about the implementation ? Thanks.

shankar1729 commented 2 years ago

See Phys. Rev. B 87, 165122 (2013) or arXiv:1302.6204. Let me know if your question is about the concept of the implementation or meaning of parameters.

MSJavaScript commented 2 years ago

See Phys. Rev. B 87, 165122 (2013) or arXiv:1302.6204. Let me know if your question is about the concept of the implementation or meaning of parameters.

Thanks, I'll read it carefully. I also read some papers about Coulomb cutoff, for example, Rozzi, et,al. Phys. Rev. B 73, 205119 (2006). But I am confused about how to deal with the singularity of truncated Coulomb kernal. AAA

As you can see, when cutoff distance tends to infinity, the kernal can't tend to 4pi/G^2.

shankar1729 commented 2 years ago

When the cutoff distance tends to infinity, the kernel will tend to 4 pi/G^2 for all G except G = 0. Here's the simplest way to think about it: all of these cutoff kernels integrate 1/r with a Fourier phase factor over some subset of all space. When you take that subset tending to the full space, the results should match, and this will happen for all finite G, where the integral over all space is well-defined.

However for G=0, the integral of 1/r over all space diverges to infinity. With the periodic Coulomb kernel, this is dealt with by explicitly setting the G=0 component to zero and making sure you have neutral unit cells etc. With the truncated kernels, there is no problem at finite cutoffs, and that is the advantage of these methods. However, taking the inifinite cutoff limit will reintroduce the periodic kernel problem.

Next, for the specific example you show for the slab kernel from my old Phys Rev paper, that expression already performed some simplifications compared to the Rozzi one. Specifically, that Gz is a multiple of a 2pi/L (or pi/R in Rozzi's notation) in the plane-wave expansion of charge densities expanded on a unit cell of that size, and so sin(GzR) = 0. This simplification is not valid in the infinite cutoff limit, and so if you wanted to take that limit you have to keep the form as in Rozzi's paper. But remember that the infinite limit will never actually be used in practice, so this is only an analytical check that we are considering.

Finally, for the G=0 limit of the slab kernel particularly, you are concerned that 4pi/G^2 -> +infinity, while the slab kernel in its infinite cutoff limit has a G=0 that tends to -infinity. This is because the G=0 component is not well-defined even for the slab kernel, and there is a choice of boundary condition being made here for definiteness. If I remember correctly, the choice is to make the potential due to a sheet charge at origin be zero at +/-L/2. Essentially, all the G=0 disagreements amount to different choices of dealing with an undetermined overall offset in the electrostatic potential.

Hope that helps with a more thorough understanding of this rather complicated problem!

Best, Shankar

MSJavaScript commented 2 years ago

Thanks for your kindly and detailed explanation, it helps a lot. In fact, I am trying to add Coulomb cutoff in VASP. However, the G vectors in VASP don't have a factor 2pi, so I have to use the formula in Rozzi's paper. If Gx, Gy and Gz ≠ 0, because the exp(-G_|| R) factor, the cutoff kernel is close to bare kernel, but when Gx=Gy=0, the kernel changes its sign at some Gz, that violates the potential completely and I get totally wrong result. (But I wonder why it works in OCTOPUS, it uses exactly the same formulas as those in Rozzi's paper) Anyway, thanks for your detailed answer, it partly solves my problem.