rhonmt.f90 4.3 KB
Newer Older
1 2
MODULE m_rhonmt
CONTAINS
3
  SUBROUTINE rhonmt(atoms,sphhar,we,ne,sym, acof,bcof,denCoeffs,ispin)
4 5 6 7 8 9 10
    !     *************************************************************
    !     subroutine sets up the coefficients of non-sphereical
    !     muffin-tin density                          c.l.fu
    !     *************************************************************
    USE m_gaunt,ONLY:gaunt1
    USE m_types
    IMPLICIT NONE
11 12 13 14 15 16 17
    TYPE(t_sym),       INTENT(IN)    :: sym
    TYPE(t_sphhar),    INTENT(IN)    :: sphhar
    TYPE(t_atoms),     INTENT(IN)    :: atoms
    TYPE(t_denCoeffs), INTENT(INOUT) :: denCoeffs

    INTEGER,           INTENT(IN)    :: ne, ispin

18 19 20
    COMPLEX, INTENT(IN) :: acof(:,0:,:)!(nobd,0:lmaxd* (lmaxd+2),natd)
    COMPLEX, INTENT(IN) :: bcof(:,0:,:)
    REAL,    INTENT(IN) :: we(:)!(nobd)
21

22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75
    !     ..
    !     .. Local Scalars ..
    COMPLEX cconst,cil,cmv,ci
    REAL coef
    INTEGER  :: jmem,l,lcond,lh,llp,llpmax,lm,lmp,lp,lphi,lplow,lplow0,lv
    INTEGER  :: mp,mv,na,natom,nb,nn,ns,nt,m
    !     ..
    !     ..

    ci = cmplx(0.0,1.0)
    !
    !Initialize private variables in gaunt module before parallel region
    !$      coef = gaunt1(0,0,0,0,0,0,atoms%lmaxd)

    DO ns = 1,sym%nsymt
       !$OMP PARALLEL DO DEFAULT(SHARED) &
       !$OMP&  PRIVATE(lv,jmem,mv,cmv,l,lm,mp,m,llpmax,nt,na,nb,lplow0)&
       !$OMP&  PRIVATE(lphi,lplow,lcond,lp,cil,lmp,llp,coef,cconst&
       !$OMP& ,natom,nn) 
       DO  lh = 1,sphhar%nlh(ns)
          lv = sphhar%llh(lh,ns)
          DO  jmem = 1,sphhar%nmem(lh,ns)
             mv = sphhar%mlh(jmem,lh,ns)
             cmv = conjg(sphhar%clnu(jmem,lh,ns))
             DO  l = 0,atoms%lmaxd
                m_loop: DO m = -l,l
                   lm = l* (l+1) + m
                   mp = m - mv
                   !     -----> set up the lower and upper limit of lp in such a way that
                   !     -----> lp+l+lv is even, lp<=l, and (lp,l,lv) satisfies the
                   !     -----> triangular relation
                   lplow0 = iabs(l-lv)
                   lphi = l - mod(lv,2)
                   lplow = max(lplow0,iabs(mp))
                   lcond = iabs(lphi-lplow)
                   lplow = lplow + mod(lcond,2)
                   IF (lplow.GT.lphi) CYCLE m_loop
                   DO  lp = lplow,lphi,2
                      cil = ci** (l-lp)
                      lmp = lp* (lp+1) + mp
                      IF (lmp.GT.lm) CYCLE m_loop
                      llp = (l* (l+1))/2 + lp
                      !     -----> gaunt's coefficient
                      coef = 2.*gaunt1(l,lv,lp,m,mv,mp,atoms%lmaxd)
                      IF (lmp.EQ.lm) coef = coef/2.
                      cconst = coef* (cil*cmv)
                      natom = 0
                      DO  nn = 1,atoms%ntype
                         llpmax = (atoms%lmax(nn)* (atoms%lmax(nn)+3))/2
                         IF (llp.LE.llpmax) THEN 
                            nt = natom
                            DO  na = 1,atoms%neq(nn)
                               nt = nt + 1
                               IF (atoms%ntypsy(nt).EQ.ns) THEN
76 77
                                  DO nb = 1,ne
                                     denCoeffs%uunmt(llp,lh,nn,ispin) = denCoeffs%uunmt(llp,lh,nn,ispin)&
78
                                          +we(nb)*real(cconst*acof(nb,lm,nt)*conjg(acof(nb,lmp,nt)))
79
                                     denCoeffs%ddnmt(llp,lh,nn,ispin) = denCoeffs%ddnmt(llp,lh,nn,ispin) +&
80
                                          we(nb)*real(cconst*bcof(nb,lm,nt)*conjg(bcof(nb,lmp,nt)))
81
                                     denCoeffs%udnmt(llp,lh,nn,ispin) = denCoeffs%udnmt(llp,lh,nn,ispin) +&
82
                                          we(nb)*real(cconst*acof(nb,lm,nt)*conjg(bcof(nb,lmp,nt)))
83
                                     denCoeffs%dunmt(llp,lh,nn,ispin) = denCoeffs%dunmt(llp,lh,nn,ispin) +&
84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100
                                          we(nb)*real(cconst*bcof(nb,lm,nt)*conjg(acof(nb,lmp,nt)))
                                  ENDDO
                               ENDIF
                            ENDDO
                         ENDIF
                         natom = natom + atoms%neq(nn)
                      ENDDO
                   ENDDO
                ENDDO m_loop
             ENDDO
          ENDDO
       ENDDO
       !$OMP END PARALLEL DO
    ENDDO

  END SUBROUTINE rhonmt
END MODULE m_rhonmt