Spurious dependence of eps on OMEGAMAX

[imazin2]

Here is a set of test calculations for V metals, calculating the optical dielectric function. The "test#" files are imaginary parts extracted from OUTCARs. In INCAR.test# are the corresponding INCARs. Changing parameters like SIGMA, NBANDS, NEDOS etc does not affect the result. However changing OMEGAMAX alters it dramatically, especially in the high-energy region. I emphasize that this is the IMAGINARY part, before any Kramers-Kronig, so it should be completely oblivious to OMEGAMAX...

Clipboard01.png

[ahampel]

Hi,

thank you for reaching out to us on the official VASP forum. Sorry for the late response. I had to understand first in detail what the code is doing, as this was first also a bit puzzling to me. The behavior of the code might be a bit unexpected but there seems to be no bug per-se, even could discuss whether the behavior is confusing.

What happens in the code is that first the imaginary part of epsilon is computed with a gaussian broadening (ISMEAR + SIGMA INCAR flags). Afterwards, the real part of epsilon is calculated via KK and a Lorentzian broadening (width determined by CSHIFT INCAR flag) is applied. But, at the end of this routine the imaginary part of epsilon is also replaced by this broadened KK transformed epsilon! Hence, the Gaussian and Lorentzian broadening are both applied to the imaginary part as well. The only times this is not the case is when CSHIFT is very small compared to the energy spacing of the epsilon grid. Now, comes maybe the part leads to the most confusion: OMEGAMAX does not cut off the transitions entering the calculation of epsilon but rather cuts off epsilon itself! This is different to what is done in the BSE routines of VASP. Using now this Im epsilon that is just cut off at certain OMEGAMAX and apply a Lorentzian broadening and you will suppress the signal in the way you see in your plot. I can demonstrate this with your system. I used SIGMA=0.1 and a rough kpoint mesh of 5x5x5 to only have a few delta peaks broadened. The I vary CSHIFT for two sets of calculation OMEGAMAX=20 and OMEGAMAX=3 :

cshift_scan_1e-01.png

cshift_scan_1e-02.png

cshift_scan_1e-03.png

cshift_scan_1e-04.png

As you can see lowering the lorentzian width recovers the omegamax=20 results.This shows clearly that the effect is not due to cut-off transitions, but only due to the convolution with the Lorentzian. In your calculations you used the default cshift=0.1 which is maybe a bit large here. So such large cshift will broaden the transition delta peaks beyond the omegamax and one looses their spectral weight.

I hope this makes sense, let me know if you have questions or if you see some other sources of problems. I will try to adapt the wiki to clarify things a bit.

Best,
Alex

[imazin2]

Many thanks! It does make total sense.
I thought it could have been CSHIFT, but VASP manual says quite clearly that CSHIFT applies only to the real part, so I abandoned this idea.

If I may suggest, modifying the manual to reflect what you just explained would go a long way.

Vielen Dank!
Igor

[ahampel]

Hi Igor,

yes I went ahead and changed parts of wiki/index.php/LOPTICS#Spectral_broadening and parts of wiki/index.php/OMEGAMAX to hopefully make this clearer.

Best regards,
Alex

[imazin2]

Hi Alex,
Since you've already spent time digging into optics here, maybe I could suggest a few improvements that shouldn't be too time consuming to implement but may be very helpful to the community:

1. Currently the number of frequencies is set to be 1/10th of NEDOS (which, incidentally, is not documented). It would be much better to decouple these two unrelated variables and introduce something like NOMEGA
2. Currently, there is no way to do optics without recalculating the band structure. The former is much less time consuming than the latter. Given that all information needed for doing optics (EIGENVAL, WAVECAR, WAVEDER) is already saved on the disk, it would be very helpful to be able to redo the optics with different optical parameters (OMEGAMAX, CSHIFT) without redoing the electronic structure. For instance, making LOPTICS take three options: 0 (LOPTICS=F), 1 (LOPTICS=T) and 2 (recalculate optics using existing files).
3. Currently, the code does not check optical parameters for consistency until the bands are finished, so one can wait for an hour just to get a message that ISMEAR = -5 is not compatible with OPTICS.
4. Last but not least. I have a universal optical tetrahedron integration routine (which is used, for instance, in Wien2k), which I wrote when I was a postdoc. It is completely universal - it takes two sets of four energies at the corners of a tetrahedron and calculates \int{M_k delta(E_1k-E_2k-w) [theta(E2k)-theta(E1k)].
It is lightning-fast, because the entire w array is computed simultaneously and converges much faster with k-points. If you would be willing to implement it in VASP I'd be glad to share this routine.
You can write me directly at imazin2@gmu.edu
Cheers,
Igor

[ahampel]

Hi Igor,

thank you for your feedback. We will add this as an issue in our internal issue tracker and try to improve the current optics routines. We will definitely try to fix this in the near future. Regarding point 4 I will contact you directly.

Best,
Alex

[imazin2]

Hi Alex, I received your direct mail and tried to answer twice. Did you receive it?