My Community

Hello,

I am trying to calculate a frequency-dependent dielectric tensor, and as per the instructions on the VASP wiki I am starting with a static calculation using LEPSILON=True, and then trying to converge a frequency-dependent calculation using LOPTICS=True by increasing NBANDS until the zero-frequency value equals the static values with LEPSILON=True.

I don't think I have enough bands yet because the diagonal values are at omega=0 are still quite a bit off from the LEPSILON results, (even though I have tripled the number of bands with respect to the default). But an even more concerning issue is that the offdiagonal components for the tensor with LEPSILON=true were antisymmetric, that is eps_xy=-eps_yx. However, separately from being a different magnitude, I see that with LOPTICS=True, the off-diagonal components at all frequencies include omega=0 are symmetric, that is eps_xy=eps_yx. This seems to me like it but be do to issues with the symmetry routine messing up in one of the two cases, as I can't imagine swapping the symmetry-allowed qualitative form of the tensor would simply be a convergence issue.

I attach the INCARs for both the static and frequency-dependent INCARs as well as the other input files (the calculation for LOPTICS=True starts from a converged scf WAVECAR as instructed). If you see any issues (of if I am missing something and it's not surprising that the symmetry of the off-diagonal elements switches between the two methods) I would be grateful. Let me know if you need the OUTCAR as well.

Hello,

I have an update/partial explanation of this. I think the problem might be on my end, or rather, on the pymatgen function I am using to extract the frequency-dependent matrix from the OUTCAR file. https://github.com/materialsproject/pym ... 2558-L2623

I noticed that, anyway as far as I can tell, the OUTCAR with LOPTICS only prints xy, yz and zx off-diagonal components and does not explicitly give the transport yx, zy, and xz components. Can you confirm that that's the case? From the source code of the pymatgen function read_freq_dielectric, I think this code just assumes that eps_yx =eps_xy (which is often true, but not in my case, where eps_yx=-eps_xy, as it's correctly given in the zero-frequency LEPSILON=True calculation I did.

If my interpretation is correct, how are you supposed to infer the transpose components of the tensor for LOPTICS=.TRUE. in general? Is it that you compare your results to the LEPSILON=True case, or is the full 3x3 matrix printed in one of the other output files?

Thanks a lot in advance. More general related to convergence with respect to LOPTICS=TRUE, are there any other parameters I should play around with other than NBANDS? I am surpised that after nearly quadrupling the number of bands with respect to default, I'm still not getting agreement in the zero frequency limit with LEPSILON=True quantitative values.

I consulted two colleagues. We had to try out a few things and ended up looking into literature.
Concerning the different values between LEPSILON=.TRUE. and LOPTIC=.TRUE.:
In the LOPTIC calculations one gets the real part of the dielectric function by making a Kramers-Kronig transformation of the imaginary dielectric function. To converge this transformation one needs many bands, so you could still not be converged with the number of bands, especially, because your structure has a very small band gap and hence high values for the dielectric constant. This is also a problem for the k-mesh since the structure is almost a metal and the number of k-points could also be not converged. So please try converging both the number of k-points and the number of bands (mind: it could also be that LEPSILON and LOPTIC converge at a different rate).

Concerning the symmetry and antisymmetry problem:
Non-magnetic calculations are indeed symmetric as can be found here.
In your calculation you are performing a non-collinear calculation for a magnetic compound. In that case one should expect an antisymmetric solution. For an explanation please have a look at this publication.

A snapshot from the publication that explains the behaviour is attached to this post here:

Hello,

Thank you, this is helpful. based on what you said, I think that probably the parameter I need to converge better is the in-plane (x-y) k mesh. I am at 576 bands now (the default was 108) with a 13x13x4 k-mesh, and the zz component of epsilon at omega=0 for LOPTICS=TRUE is almost exactly matching with LEPSILON=True, whereas the epsilon xx yy and xy components are still way off. But just to make sure I understand, is it ok/ correct to compare LOPTICS resuls with a denser in-plane k-grid to the LEPSILON=TRUE values with a lower k-grid? This was my understanding of what you said, i.e. the k-mesh for LOPTICS might have to be larger than for LEPSILON. Or should I increase the k-mesh for both types of calculations so that I'm always comparing with the same k-mesh?

(Of course, I understand that I ought to first verify that the values for LEPSILON=True do not change if I increase the k-mesh, but assuming they don't my question would then be whether I ought to still increase the mesh for LOPTICS).

For convergence I had one more question regarding choosing the value of CSHIFT; I understand this isn't really a convergence parameter, but it strongly effects the qualitative look, and quantitative values, of the epsilon (omega) curves; I attache pictures with both CSHIFT=0.1 and CSHIFT=0.4 the omega=0 values in both cases for in-plane dielectric components are way off the values from LEPSILON=True (see attached screenshot for this tensor), but the plot looks much nicer with CSHIFT=0.4. From the wiki though, it seems like this is a higher value than recommended.

Since CSHIFT strongly affects the values of epsilon(0) and I am in the end trying to converge these with respect to the result for LEPSILON=TRUE., can you clarify how CSHIFT ought to be chosen?

Both LOPTICS and LEPSILON calculations need to be converged with respect to k-points then you can compare them. Both methods are not guaranteed to converge at the same rate, so comparing uncoverged values will not work.