Unable to open bandgap of FeO (GGA+U & HSE06)

[sukhito teh]

Dear Staffs and Users,
I noticed that both GGA+U and HSE06 failed to open the bandgap of FeO. FYI, calculations with VASP can have similar issue, but those can usually be solved by turning off the K-Points symmetry. However,  even if I am setting use_symmetries  = False in AlgorithmParameters, the calculation still assume time reversal symmetry, so I am not sure if similar approach would solve the problem. May I know if there is any suggestion to solve this problem?
Thank you for your time.

best regards,
Sukhito

[Anders Blom]

The keyword you are looking for is force_timereversal=False (https://docs.quantumatk.com/manual/Types/MonkhorstPackGrid/MonkhorstPackGrid.html)

[sukhito teh]

Thank you for your reply. I tried rerun with this setting, but it doesn't seem to work.

[Anders Blom]

Define "doesn't work". The calculation ran fine and I don't any density of states being calculated so how can you tell there is no band gap?

You might want to try the Dielectric Dependent Hybrid model as discussed in https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8007096/ where FeO is an included example. This model is available QuantumATK using

exchange_correlation = HybridGGA.HSE06DDH

(cf. https://docs.quantumatk.com/manual/DFTLCAO.html#dielectric-dependent-hybrid-functionals)

[sukhito teh]

Define "doesn't work". The calculation ran fine and I don't any density of states being calculated so how can you tell there is no band gap?

You might want to try the Dielectric Dependent Hybrid model as discussed in https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8007096/ where FeO is an included example. This model is available QuantumATK using

exchange_correlation = HybridGGA.HSE06DDH

(cf. https://docs.quantumatk.com/manual/DFTLCAO.html#dielectric-dependent-hybrid-functionals)

I am so sorry that I didn't make it clear in my previous reply, doesn't work refer to unable to turn off the time reversal symmetry. Anyway, I tried few more approaches and these are the summary:
1. Using SG15 basis + PBEU --> no Band Gap
2. USing SG15 basis + HSE06 --> no Band Gap
3. Using Dojo basis +HSE06 --> no Band Gap
4. Using Dojo basis + DDH --> no Band Gap
5. Using Dojo basis + local DDH --> no Band Gap
6. Using Dojo basis + local DDH + ADMM --> band gap open!

I am not sure what is the explanation for this, but I am really thankful for your help. The ADMM is interesting, it speeds up the calculation for more than 10 times (per step), although it takes much larger steps to converge. Have a nice day!

[Anders Blom]

ADMM is indeed a method to speed up the calculations, but in theory should not really influence the results. It might be a random effect in the numerics, but no bandgap vs a bandgap is not a small effect. I wonder if we can solve the time-reversal issue and combine if with HSE06 (no ADMM) to see a band gap there, and then just use ADMM as intended to speed it up.

Can you share a few more details about how you determine the time-reversal didn't work?

[sukhito teh]

ADMM is indeed a method to speed up the calculations, but in theory should not really influence the results. It might be a random effect in the numerics, but no bandgap vs a bandgap is not a small effect. I wonder if we can solve the time-reversal issue and combine if with HSE06 (no ADMM) to see a band gap there, and then just use ADMM as intended to speed it up.

Can you share a few more details about how you determine the time-reversal didn't work?

Yes, I believe it is numeric random effect, right now I can't reproduce the band gap opening's results (with ADMM is on). I initially ran a calculation (with ADMM on) and able to open the band gap, I later rewrite the script with some minor changes in parameters and using the previous calculation (with band gap) as initial state, this indeed reproduce band gap solution as expected. However, I had forgotten what parameters had been used in the initial run and cant reproduce the band gap result without reading the "right" initial state, so the band gap opening is probably unrelated to ADMM. May I know if there is a way to initialize different density matrix during beginning of  calculation (without using pre-converged calculation)? Or is there a way to fix selective occupation state throughout the calculation?

I also tried a dozen of calculation with and without ADMM, in some cases (e.g. MnO), I do notice significant difference (band gap differ by about 10%),  not sure if such a difference is within reasonable range.

Regarding the time-reversal symmetry, in the bulkscf.log file, the "Diagonalization solver parallelization report" states there are 260 kpoints in 8x8x8 MK grids, so I assumed the time reversal symmetry is on (please correct me if I am wrong). 

[Anders Blom]

Interesting. I observed the same issue with the k-points, I am trying to investigate that further. The object k_point_sampling itself correctly reports 512 points, but for some reason it still gets reduced (or reported as such, at least) in the calculation.

It does speed up convergence if you initialize the spins closer to the expected converged polarization, rather than setting it to a fully polarized state as you do now. So in this way you do initialize the density matrix, partly, but there is no simple way to initialize the density matrix itself beyond that. One could run a separate calculation with looser convergence settings and maybe fewer k-points to at least give it a starting guess, but note that the basis set size must be the same.

Since ADMM is an approximation it is expected it may give slightly different results. 10% sounds a bit high though but it will depend on the size of the smaller basis set, of course.

[Anders Blom]

We have been digging into the code a bit more, and apparently it's actually not possible to turn off time-reversal symmetry for unpolarized for polarized calculations. This is consistent with the physics of the problem definition, but not so well explained in the manual. You can however turn it off by running a noncollinear spin calculation, and again that is how it's supposed to work. Unfortunately all matrices are then fully 2x2 so the calculation times go up significantly, but for 4 atoms it should still be manageable.

[Anders Blom]

We have made some internal progress on this, we will report the results soon

[Troels-Markussen]

Hi,
I have been looking into this issue. It turns out to be related to the initial occupation of the orbitals. This is controlled by the parameter "filling_method" on the BasisSet. By default this 'SphericalSymmetric' meaning that all orbitals within a particular shell (e.g. all the d-orbitals) will have the same initial occupation. However, this might trap the solution in a wrong local minimum. In your case the solution comes out with equal occupation of the xy, xz, and yz orbitals for the minority spin. However, according to this paper DOI: 10.1103/PhysRevB.97.035107 the solution should rather be a single occupied orbital (xy) for the minority spin. This you can actually achieve if you initialize the orbital occupations differently by setting  "filling_method=Anisotropic" on the Fe basis set (if you are using GGA+U this can also be set in the GUI, where you specify the U-values). By doing so you end up with a solution that indeed has a band gap - see attached PDOS plot, where I show the DOS projected on the first Fe atom on spin up and down. The script for the calculation is also attached.

[sukhito teh]

Thank you for the helps!

[Troels-Markussen]

A small update: When using filling_method=Anisotropic, HSE06 also leads to solution with a band gap.

[sukhito teh]

A small update: When using filling_method=Anisotropic, HSE06 also leads to solution with a band gap.

Yes, LDA+U also leads to result with band gap.