Line data    Source code

       1             : MODULE m_rhonmt
       2             : CONTAINS
       3           0 :   SUBROUTINE rhonmt(atoms,sphhar,we,ne,sym,eigVecCoeffs,denCoeffs,ispin)
       4             :     !     *************************************************************
       5             :     !     subroutine sets up the coefficients of non-sphereical
       6             :     !     muffin-tin density                          c.l.fu
       7             :     !     *************************************************************
       8             :     USE m_gaunt,ONLY:gaunt1
       9             :     USE m_types
      10             :     use m_constants
      11             :     IMPLICIT NONE
      12             :     TYPE(t_sym),          INTENT(IN)    :: sym
      13             :     TYPE(t_sphhar),       INTENT(IN)    :: sphhar
      14             :     TYPE(t_atoms),        INTENT(IN)    :: atoms
      15             :     TYPE(t_eigVecCoeffs), INTENT(IN)    :: eigVecCoeffs
      16             :     TYPE(t_denCoeffs),    INTENT(INOUT) :: denCoeffs
      17             : 
      18             :     INTEGER,           INTENT(IN)    :: ne, ispin
      19             : 
      20             :     REAL,    INTENT(IN) :: we(:)!(nobd)
      21             : 
      22             :     !     ..
      23             :     !     .. Local Scalars ..
      24             :     COMPLEX cconst,cil,cmv
      25             :     REAL coef
      26             :     INTEGER  :: jmem,l,lcond,lh,llp,llpmax,lm,lmp,lp,lphi,lplow,lplow0,lv
      27             :     INTEGER  :: mp,mv,na,natom,nb,nn,ns,nt,m
      28             : 
      29             :     !Initialize private variables in gaunt module before parallel region
      30           0 :     !$      coef = gaunt1(0,0,0,0,0,0,atoms%lmaxd)
      31             : 
      32           0 :     DO ns = 1,sym%nsymt
      33             :        !$OMP PARALLEL DO DEFAULT(SHARED) &
      34             :        !$OMP&  PRIVATE(lv,jmem,mv,cmv,l,lm,mp,m,llpmax,nt,na,nb,lplow0)&
      35             :        !$OMP&  PRIVATE(lphi,lplow,lcond,lp,cil,lmp,llp,coef,cconst&
      36           0 :        !$OMP& ,natom,nn)
      37             :        DO  lh = 1,sphhar%nlh(ns)
      38             :           lv = sphhar%llh(lh,ns)
      39             :           DO  jmem = 1,sphhar%nmem(lh,ns)
      40             :              mv = sphhar%mlh(jmem,lh,ns)
      41             :              cmv = conjg(sphhar%clnu(jmem,lh,ns))
      42             :              DO  l = 0,atoms%lmaxd
      43             :                 m_loop: DO m = -l,l
      44             :                    lm = l* (l+1) + m
      45             :                    mp = m - mv
      46             :                    !     -----> set up the lower and upper limit of lp in such a way that
      47             :                    !     -----> lp+l+lv is even, lp<=l, and (lp,l,lv) satisfies the
      48             :                    !     -----> triangular relation
      49             :                    lplow0 = abs(l-lv)
      50             :                    lphi = l - mod(lv,2)
      51             :                    lplow = max(lplow0,abs(mp))
      52             :                    lcond = abs(lphi-lplow)
      53             :                    lplow = lplow + mod(lcond,2)
      54             :                    IF (lplow.GT.lphi) CYCLE m_loop
      55             :                    DO  lp = lplow,lphi,2
      56             :                       cil = ImagUnit** (l-lp)
      57             :                       lmp = lp* (lp+1) + mp
      58             :                       IF (lmp.GT.lm) CYCLE m_loop
      59             :                       llp = (l* (l+1))/2 + lp
      60             :                       !     -----> gaunt's coefficient
      61             :                       coef = 2.*gaunt1(l,lv,lp,m,mv,mp,atoms%lmaxd)
      62             :                       IF (lmp.EQ.lm) coef = coef/2.
      63             :                       cconst = coef* (cil*cmv)
      64             :                       DO nn = 1,atoms%ntype
      65             :                          natom = atoms%firstAtom(nn) - 1
      66             :                          llpmax = (atoms%lmax(nn)* (atoms%lmax(nn)+3))/2
      67             :                          IF (llp.LE.llpmax) THEN
      68             :                             nt = natom
      69             :                             DO na = 1,atoms%neq(nn)
      70             :                                nt = nt + 1
      71             :                                IF (sym%ntypsy(nt).EQ.ns) THEN
      72             :                                   DO nb = 1,ne
      73             :                                      denCoeffs%uunmt(llp,lh,nn,ispin) = denCoeffs%uunmt(llp,lh,nn,ispin)&
      74             :                                           +we(nb)*real(cconst*eigVecCoeffs%abcof(nb,lm,0,nt,ispin)*conjg(eigVecCoeffs%abcof(nb,lmp,0,nt,ispin)))
      75             :                                      denCoeffs%ddnmt(llp,lh,nn,ispin) = denCoeffs%ddnmt(llp,lh,nn,ispin) +&
      76             :                                           we(nb)*real(cconst*eigVecCoeffs%abcof(nb,lm,1,nt,ispin)*conjg(eigVecCoeffs%abcof(nb,lmp,1,nt,ispin)))
      77             :                                      denCoeffs%udnmt(llp,lh,nn,ispin) = denCoeffs%udnmt(llp,lh,nn,ispin) +&
      78             :                                           we(nb)*real(cconst*eigVecCoeffs%abcof(nb,lm,0,nt,ispin)*conjg(eigVecCoeffs%abcof(nb,lmp,1,nt,ispin)))
      79             :                                      denCoeffs%dunmt(llp,lh,nn,ispin) = denCoeffs%dunmt(llp,lh,nn,ispin) +&
      80             :                                           we(nb)*real(cconst*eigVecCoeffs%abcof(nb,lm,1,nt,ispin)*conjg(eigVecCoeffs%abcof(nb,lmp,0,nt,ispin)))
      81             :                                   ENDDO
      82             :                                ENDIF
      83             :                             ENDDO
      84             :                          ENDIF
      85             :                       ENDDO
      86             :                    ENDDO
      87             :                 ENDDO m_loop
      88             :              ENDDO
      89             :           ENDDO
      90             :        ENDDO
      91             :        !$OMP END PARALLEL DO
      92             :     ENDDO
      93             : 
      94           0 :   END SUBROUTINE rhonmt
      95             : END MODULE m_rhonmt