Issues in preliminary NSCF

[Flex]

This a question about Quantum Espresso, but I as it is specific to preliminary NSCF runs for Yambo (GW), I though I would ask it here.

I work with a system that is a bit problematic, as it is a 2D monolayer with a lot of empty space, and I need a lot of very high bands to converge GW accurately, i.e. about 3000-4000 bands and two times that with spin-orbit. I also have to achieve good convergence up to the highest bands. This is what a typical input looks like :
----------------------------------------------------------
&control
calculation = 'nscf'
wf_collect= .true.
prefix = 'MoS2',
pseudo_dir='./'
/
&system
ibrav = 4, celldm(1) = 6.0203189746E+00, celldm(3) = 6.644166E+00, nat = 3, ntyp = 2,
!relax Abinit
ecutwfc = 76,
nbnd=3000,
occupations = 'smearing', smearing = 'gaussian', degauss = 0.001,
noncolin=.true.
lspinorb=.true.
starting_magnetization(1) = 0.1d0 ,
nosym = .false.
force_symmorphic= .true.
/
&electrons
diago_full_acc = .true.
diago_thr_init = 1.0e-10
!conv_thr = 1.0d-10,
/
ATOMIC_SPECIES
Mo 95.96 Mo-sp_r.upf
S 32.06 S_r.upf
ATOMIC_POSITIONS (bohr)
Mo 1.3322676296E-15 3.4758327806E+00 8.7435638920E-20
S 3.0101594873E+00 1.7379163904E+00 2.9540721569E+00
S 3.0101594873E+00 1.7379163904E+00 -2.9540721569E+00
K_POINTS AUTOMATIC
9 9 1 0 0 0
------------------------------------------------------

I use a ¨conv_thr¨ of about 1.0d-10 in the SCF to have good convergence, and I wanted to achieve a similar precision for the NSCF (because I also use a similar cut in Abinit, i.e. ¨tolwfr 1.0d-18¨). However, this means a ¨diago_thr_init¨ of about 1.0d-16, and up to now, I only managed to get to ¨diago_thr_init=1.0d-10¨, and when I tried smaller, I got the standard crash of ¨ Error in routine cdiaghg (10503): problems computing cholesky¨ which basically means it can't converge to a treshold so small.

So, I was wondering :

1) Am I simply going too far and is a cut of 1.0d-10 sufficiently converged ? How should I check the ¨quality¨ of the bands computed ?

2) Is there a way to optimize the calculation for so many bands, for example by adding unconverged buffer bands on top ? Perhaps that would help to achieve a better cut ? (if it is needed)

Thanks in advance

[Daniele Varsano]

Dear Thierry,

In my opinion, you have already enough precision considering that big precision is not really needed for higher energy bands, it is important that they are orthonormal as they enter as a closure relation.
But I could be wrong and I strongly suggest you to post this question to the qe mailing list as there you can receive advises from qe experts.
Best

Daniele

[Flex]

I did some experimentation with the number of bands, and I asked the question on the QE mailing list.

There seems to be an intrinsic problem in the way that QE is coded that makes it very unstable for computing so much bands. You can go as far as 2000-3000 bands (ecutwfc 80) if you :
- take a bigger diago_thr_init,
- use CG diagonalization
- to go a bit farther, use an artificially high ecutwfc, like 150 (more plane-waves), but it becomes obviously very slow.

Higher than 3000 bands, it becomes quickly impossible to converge (even at 1.0d-5, which seems too big). That is big problem, because in MoS2, some very high D band seem to have a contribution (I verified that with Abinit) and I probably need 3000 bands and 6000 bands with spin-orbit (if I am not mistaken, you have twice the bands) to get them in my calculation.

This a big problem for me, and as more and more people use GW, this is probably going to be a problem for more people in the future.

So, I have 2 questions :

- Is there a way to converge with less bands in Yambo, an algorithm that doesn't need the sum rule ? I saw a variable named ¨Gtermkind¨ that seems to do that. Does is affect the screening of the sigma part of the calculation, or both ? Would it help ?

- Do you think this QE band problem is severe enough to post it to QE devs or on the bug tracker ?

[Daniele Varsano]

Dear Thierry,

- Is there a way to converge with less bands in Yambo, an algorithm that doesn't need the sum rule ? I saw a variable named ¨Gtermkind¨ that seems to do that. Does is affect the screening of the sigma part of the calculation, or both ? Would it help ?

Yes, you can use terminators:

Code: Select all

Gtermkind='bg"

and

Code: Select all

XTermKind= "bg"

this is stable and helps quite a lot in converging the Gbands and bands in screening respectively.
The price to pay is a more memory usage.

- Do you think this QE band problem is severe enough to post it to QE devs or on the bug tracker ?

I do not have any suggestion on that, usually, the qe mailing list is read also by developers, and this is not strictly a bug, but it looks more a limitation of the algorithm.

Best,
Daniele

[martinspenke]

I think you should also check the convergence with respect to the selected energy window :

Code: Select all

GTermEn= 40.81708      eV    # [GW] GW terminator energy (only for kind="BG")

Bests
Martin

[Daniele Varsano]

Dear Martin,

GTermEn does not set an energy window but it is a position of a pole needed in implementing the Bruneval-Gonze method, see PRB 78, 085125 (2008). It is important it is placed at enough high energy, you can try to move it at even higher energy, but in general you can leave it at the default energy.

Best,

Daniele

[martinspenke]

Dear Daniele,

Thanks for your response.

I remember to have checked it in Abinit. The default was 2 Ha; and it was important for my case to set it to higher values.

According to Abinit page :

The keyword "gwencomp" in Abinit : sets the energy parameter used in the extrapolar approximation used to improve completeness and make the convergence against the number of bands much faster.

Best wishes
Martin

[Flex]

Dear Thierry,

- Is there a way to converge with less bands in Yambo, an algorithm that doesn't need the sum rule ? I saw a variable named ¨Gtermkind¨ that seems to do that. Does is affect the screening of the sigma part of the calculation, or both ? Would it help ?

Yes, you can use terminators:

Code: Select all

Gtermkind='bg"

and

Code: Select all

XTermKind= "bg"

this is stable and helps quite a lot in converging the Gbands and bands in screening respectively.
The price to pay is a more memory usage.

Daniele

In theory, using a terminator means using more memory. But after some testing with terminators, it seems I can converge similarly with a lot of bands (3000) and no terminator, or with a terminator (with the terminator energy converged) but with far less bands, as little as 1000 in some cases. This seems to be a win-win, because even if the terminator uses a little more memory, as I need far less bands, I actually use less memory overall, which is quite crucial in this kind of calculation.

This seems a bit too good to be true, so I will have to test it a bit more thoroughly, but it is also possibly specific to this material. I'm not sure, does this seem unrealistic ?

[Daniele Varsano]

Dear Thierry,

This is the idea of using the terminator.
Glad to know it is working properly!

Best,

Daniele