[Pw_forum] BUG in the LDA+U part of pwscf (Lixin He,Rutgers)

[Stefano de Gironcoli]

Dear All,

A bug in the LDA+U part of the PWSCF code has been reported by Lixin He
in Vanderbilt's group at Rutgers university.

The symmetrization of the occupation matrix used in LDA+U was done
incorrectly using the d2 matrices instead of their transpose.
In subroutine new_ns.F (or new_ns.f90 in the f90 version) the
symmetrization step should be performed with the following bit of fortran

                              nsnew(na,is,m1,m2) = nsnew(na,is,m1,m2) +
     +                              d2(m0,m1,isym) * nr(nb,is,m0,m00) *
     +                              d2(m00,m2,isym) / nsym


instead of the original one

                              nsnew(na,is,m1,m2) = nsnew(na,is,m1,m2) +
     +                              d2(m1,m0,isym) * nr(nb,is,m0,m00) *
     +                              d2(m2,m00,isym) / nsym

i.e. using d2^{T} instead of d2, more precisely:
change d2(m1,m0,isym)--> d2(m0,m1,isym)
       d2(m2,m00,isym)-->d2(m00,m2,isym).

Please correct your version accordingly if you plan to use LDA+U.

The bug is present in all versions containing LDA+U stuff, including
(sob) the newly released 1.2.0 version. It has been fixed in the version
under development.

A similar error was also present in the PP/projwave.f90 file where
atomic-state projected densities of states are computed (and symmetrized).

The original code was all right for crystals that only have mirror
symmetries or for which symmetries do not exchanges equivalent atoms
in the unit cell. But it could fail for crystals that have rotational
symmetries that exchange atoms.

Thanks to Lixin He for reporting the problem (and the solution !!!)

Stefano de Gironcoli