tlo.f90 10.2 KB
 Markus Betzinger committed Apr 26, 2016 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 ``````MODULE m_tlo USE m_juDFT !*********************************************************************** ! sets up the extra t-matrix elements due to the local orbitals. ! only non=zero elements are calculated ! ! p.kurz jul. 1996 !*********************************************************************** CONTAINS SUBROUTINE tlo(atoms,sphhar, jsp,ntyp,enpara,lh0,input,vr,& flo,f,g,usdus,uuilon,duilon,ulouilopn, tlmplm ) ! !*************** ABBREVIATIONS ***************************************** ! tuulo : t-matrix element of the lo and the apw radial fuction ! tdulo : t-matrix element of the lo and the energy derivativ of ! the apw radial fuction ! tuloulo : t-matrix element of two los !c*********************************************************************** ! USE m_intgr, ONLY : intgr3 USE m_gaunt, ONLY: gaunt1 USE m_types IMPLICIT NONE TYPE(t_input),INTENT(IN) :: input TYPE(t_sphhar),INTENT(IN) :: sphhar TYPE(t_atoms),INTENT(IN) :: atoms TYPE(t_usdus),INTENT(IN) :: usdus TYPE(t_tlmplm),INTENT(INOUT):: tlmplm TYPE(t_enpara),INTENT(IN) :: enpara ! .. ! .. Scalar Arguments .. INTEGER, INTENT (IN) :: jsp,ntyp ,lh0 ! .. ! .. Array Arguments .. REAL, INTENT (IN) :: vr(atoms%jmtd,0:sphhar%nlhd) REAL, INTENT (IN) :: f(atoms%jmtd,2,0:atoms%lmaxd),g(atoms%jmtd,2,0:atoms%lmaxd) REAL, INTENT (IN) :: flo(atoms%jmtd,2,atoms%nlod) REAL, INTENT (IN) :: uuilon(atoms%nlod,atoms%ntypd),duilon(atoms%nlod,atoms%ntypd) REAL, INTENT (IN) :: ulouilopn(atoms%nlod,atoms%nlod,atoms%ntypd) ! .. ! .. Local Scalars .. COMPLEX ci,cil `````` Daniel Wortmann committed Apr 27, 2016 43 `````` INTEGER i,l,lh,lm ,lmin,lmp,lo,lop,loplo,lp,lpmax,lpmax0,lpmin,lpmin0,lpp ,mem,mp,mpp,m,lmx,mlo,mlolo `````` Markus Betzinger committed Apr 26, 2016 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 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 `````` ! .. ! .. Local Arrays .. REAL x(atoms%jmtd),ulovulo(atoms%nlod*(atoms%nlod+1)/2,lh0:sphhar%nlhd) REAL uvulo(atoms%nlod,0:atoms%lmaxd,lh0:sphhar%nlhd),dvulo(atoms%nlod,0:atoms%lmaxd,lh0:sphhar%nlhd) ! .. ci = CMPLX(0.0,1.0) DO lo = 1,atoms%nlo(ntyp) l = atoms%llo(lo,ntyp) DO lp = 0,atoms%lmax(ntyp) lmin = ABS(lp-l) ! lmin = lp - l lmx = lp + l DO lh = lh0,sphhar%nlh(atoms%ntypsy(atoms%nat)) lpp = sphhar%llh(lh,atoms%ntypsy(atoms%nat)) IF ((MOD(l+lp+lpp,2).EQ.1) .OR. (lpp.LT.lmin) .OR.& (lpp.GT.lmx)) THEN uvulo(lo,lp,lh) = 0.0 dvulo(lo,lp,lh) = 0.0 ELSE DO i = 1,atoms%jri(ntyp) x(i) = (f(i,1,lp)*flo(i,1,lo)+ f(i,2,lp)*flo(i,2,lo))*vr(i,lh) END DO CALL intgr3(x,atoms%rmsh(:,ntyp),atoms%dx(ntyp),atoms%jri(ntyp),uvulo(lo,lp,lh)) DO i = 1,atoms%jri(ntyp) x(i) = (g(i,1,lp)*flo(i,1,lo)+ g(i,2,lp)*flo(i,2,lo))*vr(i,lh) END DO CALL intgr3(x,atoms%rmsh(:,ntyp),atoms%dx(ntyp),atoms%jri(ntyp),dvulo(lo,lp,lh)) END IF END DO END DO END DO loplo = 0 DO lop = 1,atoms%nlo(ntyp) lp = atoms%llo(lop,ntyp) DO lo = 1,lop l = atoms%llo(lo,ntyp) loplo = loplo + 1 IF (loplo>size(ulovulo,1)) CALL juDFT_error("loplo too large!!!" ,calledby ="tlo") DO lh = lh0,sphhar%nlh(atoms%ntypsy(atoms%nat)) lpp = sphhar%llh(lh,atoms%ntypsy(atoms%nat)) lmin = ABS(lp - l) lmx = lp + l IF ((MOD(l+lp+lpp,2).EQ.1).OR.(lpp.LT.lmin).OR.(lpp.GT.lmx)) THEN ulovulo(loplo,lh) = 0.0 ELSE DO i = 1,atoms%jri(ntyp) x(i) = (flo(i,1,lop)*flo(i,1,lo)+flo(i,2,lop)*flo(i,2,lo))*vr(i,lh) END DO CALL intgr3(x,atoms%rmsh(:,ntyp),atoms%dx(ntyp),atoms%jri(ntyp),ulovulo(loplo,lh)) END IF END DO END DO END DO !---> generate the different t matrices !---> but first initialize them ( done in eigen ) ! !---> generate the t-matrices. 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 !---> loop over the local orbitals `````` Daniel Wortmann committed Apr 27, 2016 107 `````` mlo=sum(atoms%nlo(:ntyp-1)) `````` Markus Betzinger committed Apr 26, 2016 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 `````` DO lo = 1,atoms%nlo(ntyp) l = atoms%llo(lo,ntyp) DO m = -l,l !---> loop over the lattice harmonics DO lh = lh0,sphhar%nlh(atoms%ntypsy(atoms%nat)) lpp = sphhar%llh(lh,atoms%ntypsy(atoms%nat)) 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 mem = 1,sphhar%nmem(lh,atoms%ntypsy(atoms%nat)) mpp = sphhar%mlh(mem,lh,atoms%ntypsy(atoms%nat)) 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 cil = ((ci** (l-lp))*sphhar%clnu(mem,lh,atoms%ntypsy(atoms%nat)))* gaunt1(lp,lpp,l,mp,mpp,m,atoms%lmaxd) `````` Daniel Wortmann committed Apr 27, 2016 131 132 133 134 `````` tlmplm%tuulo(lmp,m,lo+mlo,jsp) = & tlmplm%tuulo(lmp,m,lo+mlo,jsp) + cil*uvulo(lo,lp,lh) tlmplm%tdulo(lmp,m,lo+mlo,jsp) = & tlmplm%tdulo(lmp,m,lo+mlo,jsp) + cil*dvulo(lo,lp,lh) `````` Markus Betzinger committed Apr 26, 2016 135 136 137 138 139 140 141 `````` END DO END DO END DO END DO END DO !---> generate the t-matrix including two local orbitals for lo' >= lo !---> loop over lo' `````` Daniel Wortmann committed Apr 27, 2016 142 `````` mlolo=dot_product(atoms%nlo(:ntyp-1),atoms%nlo(:ntyp-1)+1)/2 `````` Markus Betzinger committed Apr 26, 2016 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 `````` DO lop = 1,atoms%nlo(ntyp) lp = atoms%llo(lop,ntyp) DO mp = -lp,lp !---> loop over the lattice harmonics DO lh = lh0,sphhar%nlh(atoms%ntypsy(atoms%nat)) lpp = sphhar%llh(lh,atoms%ntypsy(atoms%nat)) DO mem = 1,sphhar%nmem(lh,atoms%ntypsy(atoms%nat)) mpp = sphhar%mlh(mem,lh,atoms%ntypsy(atoms%nat)) m = mp - mpp !---> loop over lo DO lo = 1,lop l = atoms%llo(lo,ntyp) loplo = ((lop-1)*lop)/2 + lo IF ((ABS(l-lpp).LE.lp) .AND. (lp.LE. (l+lpp)) .AND.& (MOD(l+lp+lpp,2).EQ.0) .AND. (ABS(m).LE.l)) THEN cil = ((ci** (l-lp))*sphhar%clnu(mem,lh,atoms%ntypsy(atoms%nat)))* gaunt1(lp,lpp,l,mp,mpp,m,atoms%lmaxd) `````` Daniel Wortmann committed Apr 27, 2016 159 `````` tlmplm%tuloulo(mp,m,loplo+mlolo,jsp) = tlmplm%tuloulo(mp,m,loplo+mlolo,jsp) + cil*ulovulo(loplo,lh) `````` Markus Betzinger committed Apr 26, 2016 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 `````` END IF END DO END DO END DO END DO END DO !---> add the diagonal terms from the muffin-tin hamiltonian. these !---> terms have to be made hermitian. if second variation is switched !---> on, the t-matrices contain only the contributions from the !---> non-spherical hamiltonian. IF (.NOT.input%secvar) THEN DO lo = 1,atoms%nlo(ntyp) l = atoms%llo(lo,ntyp) DO m = -l,l lm = l* (l+1) + m `````` Daniel Wortmann committed Apr 27, 2016 175 `````` tlmplm%tuulo(lm,m,lo+mlo,jsp) = tlmplm%tuulo(lm,m,lo+mlo,jsp) + 0.5 * usdus%uulon(lo,ntyp,jsp) *& `````` Markus Betzinger committed Apr 26, 2016 176 `````` ( enpara%el0(l,ntyp,jsp)+enpara%ello0(lo,ntyp,jsp) ) `````` Daniel Wortmann committed Apr 27, 2016 177 `````` tlmplm%tdulo(lm,m,lo+mlo,jsp) = tlmplm%tdulo(lm,m,lo+mlo,jsp) + 0.5 * usdus%dulon(lo,ntyp,jsp) *& `````` Markus Betzinger committed Apr 26, 2016 178 179 `````` ( enpara%el0(l,ntyp,jsp)+enpara%ello0(lo,ntyp,jsp) ) + 0.5 * usdus%uulon(lo,ntyp,jsp) IF (atoms%ulo_der(lo,ntyp).GE.1) THEN `````` Daniel Wortmann committed Apr 27, 2016 180 181 `````` tlmplm%tuulo(lm,m,lo+mlo,jsp) = tlmplm%tuulo(lm,m,lo+mlo,jsp) + 0.5 * uuilon(lo,ntyp) tlmplm%tdulo(lm,m,lo+mlo,jsp) = tlmplm%tdulo(lm,m,lo+mlo,jsp) + 0.5 * duilon(lo,ntyp) `````` Markus Betzinger committed Apr 26, 2016 182 183 184 `````` ENDIF !+apw_lo IF (atoms%l_dulo(lo,ntyp)) THEN `````` Daniel Wortmann committed Apr 27, 2016 185 186 `````` tlmplm%tuulo(lm,m,lo+mlo,jsp) = tlmplm%tuulo(lm,m,lo+mlo,jsp) + 0.5 tlmplm%tdulo(lm,m,lo+mlo,jsp) = 0.0 `````` Markus Betzinger committed Apr 26, 2016 187 188 189 190 191 192 193 194 195 `````` ENDIF !+apw_lo END DO END DO DO lop = 1,atoms%nlo(ntyp) lp = atoms%llo(lop,ntyp) DO lo = atoms%lo1l(lp,ntyp),lop loplo = ((lop-1)*lop)/2 + lo DO m = -lp,lp `````` Daniel Wortmann committed Apr 27, 2016 196 `````` tlmplm%tuloulo(m,m,loplo+mlolo,jsp) = tlmplm%tuloulo(m,m,loplo+mlolo,jsp) + 0.5* (enpara%ello0(lop,ntyp,jsp)+& `````` Markus Betzinger committed Apr 26, 2016 197 198 199 200 201 202 203 204 205 `````` enpara%ello0(lo,ntyp,jsp))* usdus%uloulopn(lop,lo,ntyp,jsp) + 0.5* (ulouilopn(lop,lo,ntyp) +& ulouilopn(lo,lop,ntyp)) END DO END DO END DO END IF END SUBROUTINE tlo END MODULE m_tlo``````