rhonmtlo.f90 6.78 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

7 8 9 10 11 12 13 14 15 16 17
MODULE m_rhonmtlo
  !
  !***********************************************************************
  ! This subroutine is the equivalent of rhomt for the local orbital
  ! contributions to the charge.
  ! acnmt, bcnmt, ccnmt are the equivalents of uunmt, ddnmt, udnmt dunmt
  ! in rhonmt
  ! p.kurz sept. 1996
  !***********************************************************************
  !
CONTAINS
18
  SUBROUTINE rhonmtlo(atoms,sphhar,ne,we,eigVecCoeffs,denCoeffs,ispin)
19 20
    USE m_gaunt,ONLY:gaunt1
    USE m_types
21
    use m_constants
22

23
    IMPLICIT NONE
24

25 26 27 28
    TYPE(t_sphhar),       INTENT(IN)    :: sphhar
    TYPE(t_atoms),        INTENT(IN)    :: atoms
    TYPE(t_eigVecCoeffs), INTENT(IN)    :: eigVecCoeffs
    TYPE(t_denCoeffs),    INTENT(INOUT) :: denCoeffs
29 30 31

    INTEGER, INTENT (IN) :: ne, ispin

32
    REAL,    INTENT (IN) :: we(:)!(nobd)
33

34
    !     .. Local Scalars ..
35
    COMPLEX cmv,fact,cf1
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 76 77 78
    INTEGER i,jmem,l,lh,lmp,lo,lop,lp,lpmax,lpmax0,lpmin,lpmin0,m,lpp ,mp,mpp,na,neqat0,nn,ntyp
    !     ..
    !     ..

    !---> for optimal performance consider only
    !---> those combinations of l,l',l'',m,m',m'' that satisfy the three
    !---> conditions for non-zero gaunt-coeff. i.e.
    !---> |l - l''| <= l' <= l + l'' (triangular condition)
    !---> m' + m'' = m and l + l' + l'' even

    neqat0 = 0
    DO ntyp = 1,atoms%ntype
       !--->    loop over the lattice harmonics
       DO lh = 1,sphhar%nlh(atoms%ntypsy(neqat0+1))
          lpp = sphhar%llh(lh,atoms%ntypsy(neqat0+1))
          DO jmem = 1,sphhar%nmem(lh,atoms%ntypsy(neqat0+1))
             mpp = sphhar%mlh(jmem,lh,atoms%ntypsy(neqat0+1))
             cmv = CONJG(sphhar%clnu(jmem,lh,atoms%ntypsy(neqat0+1)))
             DO lo = 1,atoms%nlo(ntyp)
                l = atoms%llo(lo,ntyp)
                lpmin0 = ABS(l-lpp)
                lpmax0 = l + lpp
                !--->             check that lpmax is smaller than the max l of the
                !--->             wavefunction expansion at this atom
                lpmax = MIN(lpmax0,atoms%lmax(ntyp))
                !--->             make sure that l + l'' + lpmax is even
                lpmax = lpmax - MOD(l+lpp+lpmax,2)
                DO m = -l,l

                   !--->                add flapw - local orbital cross-terms

                   !--->                add terms containing gaunt1(l,lp,lpp,m,mp,mpp)
                   !--->                note that gaunt1(l,lp,lpp,m,mp,mpp) computes the
                   !--->                integral of conjg(y(l,m))*y(lp,mp)*y(lpp,mpp),
                   !--->                however, since the gaunt coef. are real, this is
                   !--->                the same as int. y(l,m)*conjg(y(lp,mp)*y(lpp,mpp))
                   mp = m - mpp
                   lpmin = MAX(lpmin0,ABS(mp))
                   !--->                make sure that l + l'' + lpmin is even
                   lpmin = lpmin + MOD(ABS(lpmax-lpmin),2)
                   !--->                loop over l'
                   DO lp = lpmin,lpmax,2
                      lmp = lp* (lp+1) + mp
79
                      fact = cmv* (ImagUnit** (l-lp))*gaunt1(l,lp,lpp,m,mp,mpp,atoms%lmaxd)
80 81 82 83
                      na = neqat0
                      DO nn = 1,atoms%neq(ntyp)
                         na = na + 1
                         DO i = 1,ne
84 85 86 87 88
                            cf1 = fact *  eigVecCoeffs%ccof(m,i,lo,na,ispin)
                            denCoeffs%acnmt(lp,lo,lh,ntyp,ispin) = denCoeffs%acnmt(lp,lo,lh,ntyp,ispin) +&
                                                                   we(i) * REAL(cf1 * CONJG(eigVecCoeffs%acof(i,lmp,na,ispin)) )
                            denCoeffs%bcnmt(lp,lo,lh,ntyp,ispin) = denCoeffs%bcnmt(lp,lo,lh,ntyp,ispin) +&
                                                                   we(i) * REAL(cf1 * CONJG(eigVecCoeffs%bcof(i,lmp,na,ispin)) )
89 90 91 92 93 94 95 96 97 98 99 100
                         END DO
                      END DO
                   END DO

                   !--->                add terms containing gaunt1(lp,l,lpp,mp,m,mpp)
                   mp = m + mpp
                   lpmin = MAX(lpmin0,ABS(mp))
                   !--->                make sure that l + l'' + lpmin is even
                   lpmin = lpmin + MOD(ABS(lpmax-lpmin),2)
                   !--->                loop over l'
                   DO lp = lpmin,lpmax,2
                      lmp = lp* (lp+1) + mp
101
                      fact = cmv* (ImagUnit** (lp-l))*gaunt1(lp,l,lpp,mp,m,mpp,atoms%lmaxd)
102 103 104 105
                      na = neqat0
                      DO nn = 1,atoms%neq(ntyp)
                         na = na + 1
                         DO i = 1,ne
106 107 108 109 110
                            cf1 = fact * CONJG(eigVecCoeffs%ccof(m,i,lo,na,ispin))
                            denCoeffs%acnmt(lp,lo,lh,ntyp,ispin) = denCoeffs%acnmt(lp,lo,lh,ntyp,ispin) +&
                                                                   we(i) * REAL(cf1 * eigVecCoeffs%acof(i,lmp,na,ispin) )
                            denCoeffs%bcnmt(lp,lo,lh,ntyp,ispin) = denCoeffs%bcnmt(lp,lo,lh,ntyp,ispin) +&
                                                                   we(i) * REAL(cf1 * eigVecCoeffs%bcof(i,lmp,na,ispin) )
111 112 113 114 115 116 117 118 119 120 121
                         END DO
                      END DO
                   END DO

                   !--->                add local orbital - local orbital terms
                   DO lop = 1,atoms%nlo(ntyp)
                      lp = atoms%llo(lop,ntyp)

                      !--->                   add terms containing gaunt1(l,lp,lpp,m,mp,mpp)
                      mp = m - mpp
                      IF ((ABS(l-lpp).LE.lp) .AND.(lp.LE. (l+lpp)) .AND.(MOD(l+lp+lpp,2).EQ.0) .AND.(ABS(mp).LE.lp)) THEN
122
                         fact = cmv* (ImagUnit** (l-lp))*gaunt1(l,lp,lpp,m,mp,mpp,atoms%lmaxd)
123 124 125 126
                         na = neqat0
                         DO nn = 1,atoms%neq(ntyp)
                            na = na + 1
                            DO i = 1,ne
127
                               denCoeffs%ccnmt(lop,lo,lh,ntyp,ispin) =&
128 129
                                  denCoeffs%ccnmt(lop,lo,lh,ntyp,ispin) +&
                                  we(i) * REAL(fact * CONJG(eigVecCoeffs%ccof(mp,i,lop,na,ispin))*eigVecCoeffs%ccof(m,i,lo,na,ispin))
130 131 132 133 134 135 136 137 138 139 140 141 142 143
                            END DO
                         END DO
                      END IF

                   END DO
                END DO
             END DO
          END DO
       END DO
       neqat0 = neqat0 + atoms%neq(ntyp)
    END DO

  END SUBROUTINE rhonmtlo
END MODULE m_rhonmtlo