[Wien] Q: charge lost in SO calculation depending on LAPW / APW+lo

[Eeuwe Sieds Zijlstra]

Dear all,

For bismuth (Bi) I found the following behavior, which I do not understand:

When 6s states were described by LAPW's and 6p states by APW+lo's
and when the spin-orbit (SO) interaction was included using a relativistic 
local p1/2 orbital, 0.33 electrons were lost as can be seen from the 
following output of MIXER:

:NEC01: NUCLEAR AND ELECTRONIC CHARGE   83.00000  82.67338   1.00395

This charge disappears from the p channel as can be seen when QTL is compared 
for the above case and a case where the 6s states are also described by 
APW+lo's:

:QTL001: 1.2922 1.3993 9.8302 0.0642 ....  (above case)
:QTL001: 1.2631 1.7246 9.8203 0.0643 .... (when 6s is described by APW+lo's)

No charge was lost when 6s was described by APW+lo's, when 6s and 6p were both 
described by LAPW 's (but this calculation is unstable, so I did just one 
iteration), or when no p1/2 orbital was included. 

I checked that there are no linear dependencies in the second-variational 
secular equation (the overlap matrix is positive definite) and that this 
behavior does not depend on Emax or the choice of basis for the d states, 
i.e., APW+lo or LAPW. There were no QTL-B warnings or other warnings that I 
am aware of.

When no p1/2 orbital orbital was included, the total energy was insensitive to 
whether LAPW or APW+lo was used for 6s and 6p, so in principle these choices 
form an equally good basis for LAPW1

I found this behavior for both the A7 structure (with 2 atoms per unit cell) 
and for simple cubic Bi [space group Pm-3m (number 221), a = 6.18 Bohr,
Gamma point only with temperature smearing, Eval=0.001 Ry]. The above data are 
for the latter calculation, which is simpler. For this calculation my files 
case.in1 and case.inso are:

********* case.in1 *************
WFFIL        (WFPRI, SUPWF)
  7.00       10    4 (R-MT*K-MAX; MAX L IN WF, V-NMT
 .07680   4   0      global e-param with N other choices, napw
 0   -0.364     0.010 CONT 0
 1    0.205     0.000 CONT 1
 2   -1.273     0.010 CONT 0
 2    0.113     0.000 CONT 0
K-VECTORS FROM UNIT:4   -4.0       0.9      emin/emax window
********** end of case.in1 ************

********* case.inso ***********
WFFIL
 4  4  0                      llmax,ipr,kpot
 -10.0000   0.90000           emin,emax (output energy window)
   0.  0.  1.                 direction of magnetization (lattice vectors)
  1                           number of atoms for which RLO is added
 1  0.205  0.000  CONT
 0 0 0 0 0             number of atoms for which SO is switch off;atoms
********** end of case.inso ************

Finally, I checked that this same behavior occured both with the precompiled 
WIEN2k package and when I compiled WIEN2k with the Portland Group compiler 
pgf90, with and without optimization.

My question is, whether anyone has had similar problems? Whether somebody 
could comment on its cause? I would really like to know, whether there is a 
bug in the program or whether there is some criterium according to which I 
can choose my basis set in one way but not in another way and of which I am 
not aware.

I hope that this question is not too hard. 
Any help would of course be greatly appreciated! 


Sincerely,

Eeuwe Sieds Zijlstra
Theoretical Physics
University Kassel
Germany