HSE06 calculation and band gap of Si

[CStar]

I have been using the HSE06 functional as described in the VASP manual to try to calculate the correct band gap for diamond Si. It seems that I am calculating a band gap that is too small (only ~0.65 eV), so I was wondering if anyone can tell me where I am going wrong.

These are the steps I have taken

1: DFT geometry optimization with gamma-centered 4 x 4 x 4 K-point grid

2: HSE06 electronic minimization with same KPOINT grid to generate
CHGCAR file
Here is the INCAR file for this step:

Titel
SYSTEM = Si
ISTART = 1
ICHARG = 2
!!! HF PARAMETERS
LHFCALC = .TRUE.
HFSCREEN = 0.200
ALGO = DAMPED
TIME = 0.400
PRECFOCK = Nor
NKRED = 2
AEXX = 0.25
!!!
ENCUT = 400.0
ENAUG = 800.0
LREAL = .FALSE.
NELM = 200
ISMEAR = 0
SIGMA = 0.100
IBRION = -1
NSW = 0
EDIFF = 1.0E-06
EDIFFG = -5.0E-03


3: HSE06 step with ICHARG = 11; INCAR below

Titel
SYSTEM = Si
ISTART = 1
ICHARG = 11
!!! HF PARAMETERS
LHFCALC = .TRUE.
HFSCREEN = 0.200
ALGO = DAMPED
TIME = 0.400
PRECFOCK = Nor
NKRED = 2
AEXX = 0.25
!!!
ENCUT = 400.0
ENAUG = 800.0
LREAL = .FALSE.
NELM = 200
ISMEAR = 0
SIGMA = 0.100
IBRION = -1
NSW = 0
EDIFF = 1.0E-06
EDIFFG = -5.0E-03

What I find is that the band gap seen on the DOS is not too unreasonable for the 4x4x4 grid (~1.2 eV). However, if I repeat the same steps with larger grids (6x6x6, 8x8x8, 10x10x10) then there is a clear trend that the gap is decreasing, and is only around 0.65 eV for the 10x10x10 grid. Shouldn't I be expecting a value that is close to 1.1 eV? Can anyone see where I am going wrong?

[sankh]

Well your INCAR files are faulty in many ways. Exact exchange energy is not determined through the local charges but through density matrix over KS orbitals, and hence you cannot use ICHARG = 11/1. I presume that you want to substitute 25% of short local PBE exchange by same amount of short range exact exchange, and if so, you should specify GGA = RP. For post-relaxation single point HSE06 calculations, use of NKRED is not advised, and in those calculations you could use Tetrahedron integration with Blochl correction for BZ (except for Band structure calculations) (Ismear = -5, !Sigma). If your cell has converged well, you can safely use the converged wavefunctions through setting ICHARG = 0 as starting guess.

If you experience convergence problems in HSE single point post processing calculation with Tetrahedron method (I hope not), switch back to Gaussian integration method with carefully chosen smearing width.

There have been exhaustive discussions in this forum about HSE calculations, please go through them, as well.
Cheers,
Sankh

[CStar]

Thanks sankh, your reply seems to have helped a lot. I had not realized it was not reasonable to use ICHARG=11 for the HSE calculation, but of course that makes sense. One last question though: is there any reason to specify GGA=RP when I am performing a HSE calculation? If I am using a PBE PP then would it not make sense to not set the GGA flag at all? Or is there some reason why the revised version is better for this type of calculation? Thanks again.