[Wien] Problems of HCP Tb

[Stefaan Cottenier]

>     I am running wien version wien2k 08_03

(consider upgrading -- this version is almost 2 years old, and you miss
quite some features and fixes)

>     I have been trying to simulate hcp Tb for weeks, no matter whatever 
> I did (including change the parameters and use different method such as 
> GGA+U or SO), *I still couldn't get the right magnetic moment of Tb 
> atom, and the convergence criterion was met.* The experimental magnetic 
> moment of Tb is around 10 uB, but the largest magnetic moment I got by 
> wien2k_08.3 was around 7.6 uB (GGA+U+SO).

Let us first estimate what to expect (in a free atom):

The free atom configuration of Tb is 6H_15/2 (Xe-4f^96s^2): 7 f-up and 2
f-dn electrons. That gives a spin moment of 7-2=5 mu_B. The fully
occupied up-shell does not contribute to the orbital moment. The two
dn-electron occupy the m=3 and m=2 orbitals (+3, not -3, due to Hund's
3th rule for the second half of the lanthanide series), which adds 3+2
mu_B. That gives a total (f-)moment of 5+3+2=10 mu_B per Tb-atom, which
is the experimental moment which you quote (small modifications are
possible because the solid state configuration is slightly different
from the free atom configuration).

LDA/GGA will not work, because it will put the partly occupied 7-dn
shell at E_F, in the wrong way. Hence, you are right to go to LDA+U. But
with LDA+U, you can stabilize multiple configurations of the f-shell,
and you might have to try different dmat-occupations to find the 'right'
one.

The main reason why your moment is (apparently) too small, is probably
this one:

> *************************************************************
>     *case.indmc*
> *************************************************************
> -15.                      Emin cutoff energy
>  1                       number of atoms for which density matrix is 
> calculated
>  1  1  3      index of 1st atom, number of L's, L1
>  0 0           r-index, (l,s)index

If you want to calculate the orbital moment explicitly, you have to
change the last line to '1 3' *after* convergence, and run 'x lapwdm -c
-so -up' to find the orbital moment in case.scfdmup (although for LDA+U,
there is also a :ORB line in case.scf). If you didn't do that, then for
sure your total moment is too small.

Another reason could be this one:

> *************************************************************
>     *case.inorb  (here, the value of u varies from 0.0 to 0.7 Ry)*
> *************************************************************
>   1  1  0                     nmod, natorb, ipr
> PRATT  1.0                    BROYD/PRATT, mixing
>   1 1 3                          iatom nlorb, lorb
>   1                              nsic 0..AFM, 1..SIC, 2..HFM
>    0.10 0.00        U J (Ry)   Note: we recommend to use U_eff = U-J and J=0

For lanthanides, U=0.70 Ry is a typical value. 0.10 is way too small,
and will give you almost LDA-like results, with too small moments (for
GGA+U, you might need rather 0.40 or so).

> *************************************************************
> other parameters: RMT = 2.5, Kpoints: 1000, IBZ: 76
> *************************************************************

1000 k-points is pretty small, did you test convergence with respect to
k-mesh? (8000-15000 is more plausible)

>     (2)About GGA+U, it need case.indm or case.indmc ? runsp -orb or 
> runsp -orbc ?

Use GGA+U preferably in connection with spin-orbit coupling. Then you
need case.indmc.

>        About GGA+U+SO, there are *two Tb atoms which are at the 
> unequivalent position*, how to change the case.inso files ?

? No, there should be still 2 inquivalent Tb atoms.

>     I know havey rare earth elements are hard to deal with , and I 
> really need someone who can give me some enlightenment. Any suggestion 
> will be greatly appreciated.

You can find more details in

http://dx.doi.org/10.1103/PhysRevB.74.014409

and

http://dx.doi.org/10.1103/PhysRevB.77.155101

Stefaan