2D Coulomb cutoff

[tolsen]

Hi

A have a few technical questions regarding the coulomb cutoff in 2D. I have been using this for a few systems and it seems to work very well with the random integration method applied to G=(0,0,0). I guess you use the analytical expression from Phys. Rev. B 73, 205119 (2006)?

However, I am a bit curious how it is done since the coulomb potential diverges at q=0 whenever the components of G in the periodic directions vanish. I would thus think that it is nessecary to use the RIM method for all the G vectors with vanishing components G_p in the parallel direction at q=0. On the other hand, it seems to me that the divergence can be removed if the cutoff radius R is chosen to be exactly half the non-periodic cell length, since then sin(G_n R)/G_p will vanish for all G vectors (here G_n is the component of G in the non-periodic direction). Is this a good way to circumvent these divergences or is it problematic to make this choice and ignore the inherent divergence at G_p=0. In my calculations I have not chosen the cutoff to be exactly half the cell length (but close by) so I would think the divergence is there, but it works fine. How does Yambo do it?

BR
Thomas Olsen
Post Doc
Technical University of Denmark

[Daniele Varsano]

Dear Thomas,
I understand your doubts as this part of the code it is not very well documented. I promised to myself provide more documentation soon, but finally it takes time and it is not easy to find it. If understand well, you are using a 2d cutoff using the keyword "box xy"; am I right? In this case yambo does not use the analytic expression you find in the prb paper for a sheet, but perform it numerically. In particular each coulomb cutoff component it is calculated by an analytic integration over space and a numeric integration over the Bz. You can see the detailin /src/coulomb/cutoff_box.F. Now as I presume that you have a 2d sampling of the BZ, the integration (summation) over the Brill. Zone. NEED the random integration to converge. bout analytical cutoff, the 2D of the paper has never been extensively tested (does not mean it is not correct) actually I think it has been removed from the code. anyway for the analytic 2D sheet, could be implemented this recipe : PHYSICAL REVIEW B 73, 233103 2006 that it is essentially what you are suggesting, if I well understood.

Cheers,

Daniele

[tolsen]

Dear Daniele

Thanks a lot! The method of PHYSICAL REVIEW B 73, 233103 2006 was exactly what I was thinking.

Just to be clear I have a non-periodic direction along z and use

| RandGvec= 1 RL # [RIM] Coulomb interaction RS components
| CUTGeo= "box z" # [CUT] Coulomb Cutoff geometry: box/cylinder/sphere
| % CUTBox
| 0.00000 | 0.00000 | 24.56644 | # [CUT] [au] Box sides

and I have k-point sampling in the xy plane.

You write "each coulomb cutoff component it is calculated by an analytic integration over space and a numeric integration over the Bz"

So if I understand you correctly, you take the analytical expression from Phys. Rev. B 73, 205119 (2006) with G replaced by G+q. For each q-points you then numerically integrate (average) the expression in the vicinity of that q. In this way the divergence at q=0 is also avoided because you never exactly hit q=0 in the integration? I am just wondering if this integration (at all q-points) is performed at all G vectors or just G=(0,0,0) since I have RandGvec= 1 above?

BR
Thomas Olsen
Post Doc
Technical University of Denmark

[Daniele Varsano]

Dear Thomas,

Just to be clear I have a non-periodic direction along z and use

| RandGvec= 1 RL # [RIM] Coulomb interaction RS components
| CUTGeo= "box z" # [CUT] Coulomb Cutoff geometry: box/cylinder/sphere
| % CUTBox
| 0.00000 | 0.00000 | 24.56644 | # [CUT] [au] Box sides

and I have k-point sampling in the xy plane.

correct!!

So if I understand you correctly, you take the analytical expression from Phys. Rev. B 73, 205119 (2006) with G replaced by G+q. For each q-points you then numerically integrate (average) the expression in the vicinity of that q. In this way the divergence at q=0 is also avoided because you never exactly hit q=0 in the integration? I am just wondering if this integration (at all q-points) is performed at all G vectors or just G=(0,0,0) since I have RandGvec= 1 above?

No, what Yambo do in this case is:

V^cut_G=\int_0^Zc 1/r exp(iGz)dz =
=\sum_K \int_0^Zc 1/k^2 exp(iGz)epx(iKZ)dz
I just wrote the z-component, consider you the othe 2 dimension integrated to X_c Y_c for a box, to infty for the 2D
here sum_k means \int_q \sum_K
the integration from 0 to Zc in done analitically. The summation over K (and integration over q) is done numerically.
So here consider that the quality of each component dependent on the number of G vector you sum-up (yambo consider all the Gvec you have that depend on your energy cutoff). For in the integration in BZ you need the random integration as you have a 2D sampling.
Do not know if the notation is clear, hope you catch the idea, then the calculation are simple. The suboutine is a little bit documented
and you can reconstruct it.

Cheers,

Daniele

[tolsen]

Dear Daniele

Thank a lot! It makes very good sense now.

The only thing I am still a bit puzzled about is the meaning of RandGvec since the evaluation of of V^cut_G always requires a sum over G' and integral over q' for every G? Is this because the integration is only important for certain (small) Gs?

Apologies for all these tedious questions, but I just want to make sure I understand it correctly.

BR
Thomas Olsen
Post Doc
Technical University of Denmark

[Daniele Varsano]

Dear Thomas:

The only thing I am still a bit puzzled about is the meaning of RandGvec since the evaluation of of V^cut_G always requires a sum over G' and integral over q' for every G? Is this because the integration is only important for certain (small) Gs?

You are right, anyway the problematic case is for G=0, where 3d the integral of 1/q^2 differs from the one calculated as a sum of \sum_i 1/q_i^2*VolBZ/Npoints, and for very small q it is totally overestimated (at the limit of q~0 it diverges!!). Outside the first BZ the 1/|q+G|^2 is rather smooth and the 3d integral and the 1/|q+G|*Volbz/Nq are quite the same. When performing the RIM, yambo writes in the report file
the ratio of these two integral evaluation, if at the edge of your Bz they are the same you are safe, otherwise you can include the RIM for more Gvecs,of course it depends on the shape of your Bz.
Hope now it is clear.

Cheers,

Daniele