LDA+U

LDA+U methods with different definitions of the occupation number operator [16] are available for both the collinear and non-collinear calculations by the following keyword 'scf.Hubbard.U':

    scf.Hubbard.U              on       # On|Off, default=off

It is noted that the LDA+U methods can be applied to not only LDA but also GGA. The occupation number operator is specified by the following keyword 'scf.Hubbard.Occupation':

    scf.Hubbard.Occupation     dual     # onsite|full|dual, default=dual

Among three occupation number operators, only the 'dual' operator satisfies a sum rule that the trace of occupation number matrix gives the total number of electrons which is the most primitive conserved quantity in a Hubbard model. For the details of the operator 'onsite', 'full', and 'dual', see Ref. [16]. The effective U-value in eV on each orbital of species defined by

   <Definition.of.Atomic.Species
    Ni  Ni6.0S-s2p2d2     Ni_CA11S
    O   O5.0-s2p2d1       O_CA11
   Definition.of.Atomic.Species>

is specified by

   <Hubbard.U.values                 #  eV
    Ni  1s 0.0 2s 0.0 1p 0.0 2p 0.0 1d 4.0 2d 0.0
    O   1s 0.0 2s 0.0 1p 0.0 2p 0.0 1d 0.0
   Hubbard.U.values>

The beginning of the description must be $<$Hubbard.U.values, and the last of the description must be Hubbard.U.values$>$. For all the basis orbitals, you have to give an effective U-value in eV in the above format. The '1s' and '2s' mean the first and second s-orbital, and the number behind '1s' is the effective U-value for the first s-orbital. The same rule is applied to p- and d-orbitals. As an example of the LDA+U calculation, the density of states for a nickel monoxide bulk is shown for cases with an effective U-value of 0 and 4 (eV) for d-orbitals of Ni in Fig. 27, where the input file is Crys-NiO.dat in the directory 'work'. We see that the gap increases due to the introduction of a Hubbard term on the d-orbitals. The occupation number for each orbital is output to *.out file in the same form as that of decomposed Mulliken populations which starts from the title 'Occupation Number in LDA+U' as follows:

***********************************************************
***********************************************************
       Occupation Number in LDA+U and Constraint DFT

    Eigenvalues and eigenvectors for a matrix consisting
           of occupation numbers on each site
***********************************************************
***********************************************************

    1   Ni

     spin= 0

  Sum =  8.591857905308

                    1       2       3       4       5       6       7       8
  Individual     -0.0024  0.0026  0.0026  0.0038  0.0051  0.0051  0.0888  0.0950

  s           0   0.1671  0.0005 -0.0006  0.0040  0.0000  0.0005 -0.0124  0.0000
  s           1  -0.9856 -0.0030  0.0039 -0.0227 -0.0000 -0.0072  0.0066  0.0000
  px          0   0.0010  0.0004  0.0011 -0.0131  0.0004  0.0001 -0.0261 -0.0291
  py          0   0.0010  0.0006 -0.0008 -0.0130  0.0000  0.0009 -0.0271 -0.0000
  pz          0   0.0010 -0.0012 -0.0001 -0.0131 -0.0004  0.0001 -0.0261  0.0291
  px          1   0.0067  0.0023  0.0066 -0.0792 -0.0161  0.0123  0.5594  0.7062
  py          1   0.0068  0.0041 -0.0053 -0.0801 -0.0000 -0.0162  0.5797  0.0002
  pz          1   0.0067 -0.0070 -0.0005 -0.0792  0.0161  0.0123  0.5594 -0.7063
  d3z^2-r^2   0   0.0002 -0.0781 -0.0105  0.0002  0.0023  0.0014  0.0002  0.0108
  dx^2-y^2    0   0.0004 -0.0105  0.0781  0.0004 -0.0013  0.0024  0.0003 -0.0062
  dxy         0   0.0004 -0.0009 -0.0002  0.0246 -0.0421 -0.0251  0.0794 -0.0050
  dxz         0  -0.0001  0.0008 -0.0010  0.0269  0.0000  0.0478  0.0795  0.0000
  dyz         0   0.0004  0.0004  0.0008  0.0246  0.0420 -0.0251  0.0794  0.0050
  d3z^2-r^2   1  -0.0023  0.9875  0.1327 -0.0033 -0.0262 -0.0159 -0.0001 -0.0069
  dx^2-y^2    1  -0.0040  0.1326 -0.9875 -0.0056  0.0151 -0.0275 -0.0002  0.0040
  dxy         1   0.0091  0.0233  0.0052 -0.5578  0.7055  0.4249 -0.0749  0.0157
  dxz         1   0.0189 -0.0180  0.0233 -0.5964 -0.0003 -0.7958 -0.0748 -0.0000
  dyz         1   0.0091 -0.0110 -0.0212 -0.5578 -0.7052  0.4255 -0.0749 -0.0157

                    9      10      11      12      13      14      15      16
  Individual      0.0952  0.2456  0.9902  0.9974  0.9975  1.0060  1.0060  1.0137

  s           0   0.0002  0.9859 -0.0036 -0.0001  0.0000 -0.0000  0.0000 -0.0000
  ..... 
  ...

The eigenvalues of the occupation number matrix of each atomic site correspond to the occupation number to each local state given by the eigenvector. The LDA+U functional possesses multiple minima in the degree of freedom of the orbital occupation, leading to that the SCF calculation tends to be trapped to some local minimum. To find the ground state with an orbital polarization, a way of enhancing explicitly the orbital polarization is available by the following switch :

   For collinear cases

   <Atoms.SpeciesAndCoordinates           # Unit=AU
    1  Ni   0.0       0.0       0.0       10.0  6.0  on
    2  Ni   3.94955   3.94955   0.0        6.0 10.0  on  
    3   O   3.94955   0.0       0.0        3.0  3.0  on
    4   O   3.94955   3.94955   3.94955    3.0  3.0  on
   Atoms.SpeciesAndCoordinates>

   For non-collinear cases

   <Atoms.SpeciesAndCoordinates           # Unit=AU
    1  Ni   0.0       0.0       0.0       10.0  6.0  40.0 10.0 0 0 on
    2  Ni   3.94955   3.94955   0.0        6.0 10.0  40.0 10.0 0 0 on  
    3   O   3.94955   0.0       0.0        3.0  3.0  10.0 40.0 0 0 on
    4   O   3.94955   3.94955   3.94955    3.0  3.0  10.0 40.0 0 0 on
   Atoms.SpeciesAndCoordinates>

The specification of each column can be found in the section 'Non-collinear DFT'. Since the enhancement treatment for the orbital polarization is performed on each atom, you have to set the switch for all the atoms in the specification of atomic coordinates as given above. The enhancement for the atoms switched on is applied during the first few self-consistent (SC) steps, then no more enhancement are required during the subsequent SC steps. It is also emphasized that the enhancement does not always give the ground state, and that it can work badly in some case. See Ref. [16] for the details.

Figure 27: The total density of states for up-spin in NiO bulk calculated with (a) U=0 (eV) and (b) U=4 (eV) in the LDA+U method. The input file is Crys-NiO.dat in the directory 'work'.