[gpaw-users] strange behaviour of TB09: vanishing Eg of Si

[Jens Jørgen Mortensen]

Den 17.06.2020 kl. 12.03 skrev Michal Krompiec via gpaw-users:
> Hello,
> I noticed a weird behaviour of TB09. When I do a calculation on Si with 
> default settings regarding the number of bands, I get a finite Eg as 
> expected. But when I add more bands, the band gap disappears.
> See attached code and DOS plots. What is happening here?? BTW, this 
> seems specific to TB09 (no such issue with PBE, SCAN, LDA...).

It looks like MGGA's and TB09 (which is sort of an MGGA) don't work with 
fixdensity=True.  The kinetic energy density gets initialized to zero 
and then never updated.  For now, you should avoid fixdensity=True for 
those functionals.  SCAN is an MGGA, did it really work for you?

https://gitlab.com/gpaw/gpaw/-/issues/270

Jens Jørgen

> import matplotlib as mpl
> mpl.use('Agg')
> from matplotlib import pyplot as plt
> from pathlib import Path
> 
> from ase.build import bulk
> from ase.parallel import paropen, world
> 
> from gpaw import GPAW, FermiDirac,PW, MarzariVanderbilt
> from gpaw.response.df import DielectricFunction
> 
> # Ground state calculation
> a = 5.431
> atoms = bulk('Si', 'diamond', a=a)
> 
> calc = GPAW(mode=PW(800),
>              kpts={'density': 5.0, 'gamma': True},
>              parallel={'band': 1, 'domain': 1},
>              xc='TB09',
>              occupations=MarzariVanderbilt(0.001))  # use small FD smearing
> 
> atoms.set_calculator(calc)
> atoms.get_potential_energy()  # get ground state density
> # Restart Calculation with fixed density and dense kpoint sampling
> e, dos = calc.get_dos(spin=0, npts=2001, width=0.1)
> e_f = calc.get_fermi_level()
> plt.plot(e - e_f, dos)
> plt.axis([-15, 10, None, 4])
> plt.ylabel('DOS')
> plt.savefig("dos1.png")
> calc.set(kpts={'density': 15.0, 'gamma': True},  # dense kpoint sampling
>           fixdensity=True)
> atoms.get_potential_energy()
> calc.diagonalize_full_hamiltonian(nbands=70)  # diagonalize Hamiltonian
> calc.write('si_large.gpw', 'all')  # write wavefunctions
> e, dos = calc.get_dos(spin=0, npts=2001, width=0.1)
> e_f = calc.get_fermi_level()
> plt.plot(e - e_f, dos)
> plt.axis([-15, 10, None, 4])
> plt.ylabel('DOS')
> plt.savefig("dos2.png")
> 
> Best,
> Michal Krompiec
> dos2.png
> 
> _______________________________________________
> gpaw-users mailing list
> gpaw-users at listserv.fysik.dtu.dk
> https://listserv.fysik.dtu.dk/mailman/listinfo/gpaw-users
>