abclocdn_pulay.F90 9.77 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
MODULE m_abclocdnpulay
  USE m_juDFT
CONTAINS
  SUBROUTINE abclocdn_pulay(&
       &                          atoms,sym,&
       &                          noco,ccchi,kspin,iintsp,&
       &                          con1,phase,ylm,ntyp,na,k,fgp,&
14
       &                          s,nv,ne,nbasf0,alo1,blo1,clo1,&
Daniel Wortmann's avatar
Daniel Wortmann committed
15
       &                          kvec,enough,acof,bcof,ccof,&
16
       &                          acoflo,bcoflo,aveccof,bveccof,cveccof,zMat)
17 18 19 20 21 22 23 24 25 26
    !
    !*********************************************************************
    ! for details see abclocdn; calles by to_pulay
    !*********************************************************************
    !
    USE m_types
    IMPLICIT NONE
    TYPE(t_noco),INTENT(IN)   :: noco
    TYPE(t_sym),INTENT(IN)    :: sym
    TYPE(t_atoms),INTENT(IN)  :: atoms
27
    TYPE(t_mat),INTENT(IN)    :: zMat
28 29 30 31 32 33 34 35
    !     ..
    !     .. Scalar Arguments ..
    INTEGER, INTENT (IN) :: iintsp
    INTEGER, INTENT (IN) :: k,na,ne,ntyp,nv,kspin
    REAL,    INTENT (IN) :: con1 ,s
    COMPLEX, INTENT (IN) :: phase
    !     ..
    !     .. Array Arguments ..
Daniel Wortmann's avatar
Daniel Wortmann committed
36
    INTEGER, INTENT (IN) :: nbasf0(atoms%nlod,atoms%nat) 
37
    INTEGER, INTENT (IN) :: kvec(2*(2*atoms%llod+1),atoms%nlod)
Daniel Wortmann's avatar
Daniel Wortmann committed
38 39
    REAL,    INTENT (IN) :: alo1(atoms%nlod,atoms%ntype),blo1(atoms%nlod,atoms%ntype)
    REAL,    INTENT (IN) :: clo1(atoms%nlod,atoms%ntype)
40 41
    REAL,    INTENT (IN) :: fgp(3)
    COMPLEX, INTENT (IN) :: ylm( (atoms%lmaxd+1)**2 ),ccchi(2)
Daniel Wortmann's avatar
Daniel Wortmann committed
42 43 44
    COMPLEX, INTENT (INOUT) :: acof(:,0:,:)!(nobd,0:dimension%lmd,atoms%nat)
    COMPLEX, INTENT (INOUT) :: bcof(:,0:,:)!(nobd,0:dimension%lmd,atoms%nat)
    COMPLEX, INTENT (INOUT) :: ccof(-atoms%llod:,:,:,:)!(-llod:llod,nobd,atoms%nlod,atoms%nat)
45 46
    COMPLEX, INTENT (INOUT) :: acoflo(-atoms%llod:,:,:,:)
    COMPLEX, INTENT (INOUT) :: bcoflo(-atoms%llod:,:,:,:)
Daniel Wortmann's avatar
Daniel Wortmann committed
47 48 49 50
    COMPLEX, INTENT (INOUT) :: aveccof(:,:,0:,:)!(3,nobd,0:dimension%lmd,atoms%nat)
    COMPLEX, INTENT (INOUT) :: bveccof(:,:,0:,:)!(3,nobd,0:dimension%lmd,atoms%nat)
    COMPLEX, INTENT (INOUT) :: cveccof(:,-atoms%llod:,:,:,:)!(3,-atoms%llod:llod,nobd,atoms%nlod,atoms%nat)
    LOGICAL, INTENT (OUT) :: enough(atoms%nat)
Daniel Wortmann's avatar
Daniel Wortmann committed
51
    INTEGER :: nkvec(atoms%nlod,atoms%nat)
52 53 54 55 56 57 58 59 60 61
    !     ..
    !     .. Local Scalars ..
    COMPLEX ctmp,term1
    REAL,PARAMETER:: linindq=1.0e-4,eps=1.e-30
    INTEGER i,ie,l,ll1,lm,lo ,mind,nbasf,na2,lmp,m
    LOGICAL linind
    !     ..
    !     .. Local Arrays ..
    COMPLEX clotmp(-atoms%llod:atoms%llod)
    !     ..
Daniel Wortmann's avatar
Daniel Wortmann committed
62
    nkvec=0
63 64 65 66 67 68 69 70 71 72 73
    enough(na) = .TRUE.
    term1 = con1* ((atoms%rmt(ntyp)**2)/2)*phase
    !
    !---> the whole program is in hartree units, therefore 1/wronskian is
    !---> (rmt**2)/2. the factor i**l, which usually appears in the a, b
    !---> and c coefficients, is included in the t-matrices. thus, it does
    !---> not show up in the formula above.
    !
    DO lo = 1,atoms%nlo(ntyp)
       l = atoms%llo(lo,ntyp)
       IF (.NOT.((s.LE.eps).AND.(l.GE.1))) THEN
74
          IF (sym%invsat(na).EQ.0) THEN
75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90
             IF ((nkvec(lo,na)).LT. (2*atoms%llo(lo,ntyp)+1)) THEN
                enough(na) = .FALSE.
                nkvec(lo,na) = nkvec(lo,na) + 1
                nbasf = nbasf0(lo,na) + nkvec(lo,na)
                l = atoms%llo(lo,ntyp)
                ll1 = l* (l+1)
                DO m = -l,l
                   clotmp(m) = term1*CONJG(ylm(ll1+m+1))
                END DO
                !
                IF ( kvec(nkvec(lo,na),lo) == k ) THEN
                   DO ie = 1,ne
                      DO m = -l,l
                         lm = ll1 + m
                         IF (noco%l_noco) THEN
                            IF (noco%l_ss) THEN
91
                               ctmp = clotmp(m)* ccchi(iintsp)*zMat%data_c(kspin+nbasf,ie)
92
                            ELSE
93
                               ctmp = clotmp(m)*( ccchi(1)*zMat%data_c(nbasf,ie)+ccchi(2)*zMat%data_c(kspin+nbasf,ie) )
94 95
                            ENDIF
                         ELSE
96
                            IF (zmat%l_real) THEN
97
                               ctmp = zMat%data_r(nbasf,ie)*clotmp(m)
98
                            ELSE
99
                               ctmp = zMat%data_c(nbasf,ie)*clotmp(m)
100
                            ENDIF
101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120
                         ENDIF
                         acof(ie,lm,na)     = acof(ie,lm,na) +ctmp*alo1(lo,ntyp)
                         bcof(ie,lm,na)     = bcof(ie,lm,na) +ctmp*blo1(lo,ntyp)
                         ccof(m,ie,lo,na)   = ccof(m,ie,lo,na) +ctmp*clo1(lo,ntyp)
                         acoflo(m,ie,lo,na) = acoflo(m,ie,lo,na) +ctmp*alo1(lo,ntyp)
                         bcoflo(m,ie,lo,na) = bcoflo(m,ie,lo,na) +ctmp*blo1(lo,ntyp)
                         DO i = 1,3
                            aveccof(i,ie,lm,na)=aveccof(i,ie,lm,na) +fgp(i)*ctmp*alo1(lo,ntyp)
                            bveccof(i,ie,lm,na)=bveccof(i,ie,lm,na) +fgp(i)*ctmp*blo1(lo,ntyp)
                            cveccof(i,m,ie,lo,na) =cveccof(i,m,ie,lo,na) +fgp(i)*ctmp*clo1(lo,ntyp)
                         ENDDO
                      END DO
                   END DO
                   !                    write(6,9000) nbasf,k,lo,na,
                   !     +                          (clo1(lo,ntyp)*clotmp(m),m=-l,l)
                   ! 9000               format(2i4,2i2,7(' (',e9.3,',',e9.3,')'))
                ELSE
                   nkvec(lo,na) = nkvec(lo,na) - 1
                END IF
             END IF
121
          ELSEIF (sym%invsat(na).EQ.1) THEN
122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138
             IF ((nkvec(lo,na)).LT. (2* (2*atoms%llo(lo,ntyp)+1))) THEN
                enough(na) = .FALSE.
                nkvec(lo,na) = nkvec(lo,na) + 1
                nbasf = nbasf0(lo,na) + nkvec(lo,na)
                l = atoms%llo(lo,ntyp)
                ll1 = l* (l+1)
                DO m = -l,l
                   clotmp(m) = term1*CONJG(ylm(ll1+m+1))
                END DO
                !
                IF ( kvec(nkvec(lo,na),lo) == k ) THEN
                   !                     write(*,*)'k vector nr ',k,' has been accepted'
                   DO ie = 1,ne
                      DO m = -l,l
                         lm = ll1 + m
                         IF (noco%l_noco) THEN
                            IF (noco%l_ss) THEN
139
                               ctmp = clotmp(m)*ccchi(iintsp)*zMat%data_c(kspin+nbasf,ie)
140
                            ELSE
141
                               ctmp = clotmp(m)*( ccchi(1)*zMat%data_c(nbasf,ie)+ccchi(2)*zMat%data_c(kspin+nbasf,ie) )
142 143
                            ENDIF
                         ELSE
144
                            IF (zmat%l_real) THEN
145
                               ctmp = zMat%data_r(nbasf,ie)*clotmp(m)
146
                            ELSE
147
                               ctmp = zMat%data_c(nbasf,ie)*clotmp(m)
148
                            END IF
149 150 151 152 153 154 155 156 157 158 159
                         ENDIF
                         acof(ie,lm,na) = acof(ie,lm,na) +ctmp*alo1(lo,ntyp)
                         bcof(ie,lm,na) = bcof(ie,lm,na) +ctmp*blo1(lo,ntyp)
                         ccof(m,ie,lo,na) = ccof(m,ie,lo,na) +ctmp*clo1(lo,ntyp)
                         acoflo(m,ie,lo,na) = acoflo(m,ie,lo,na) +ctmp*alo1(lo,ntyp)
                         bcoflo(m,ie,lo,na) = bcoflo(m,ie,lo,na) +ctmp*blo1(lo,ntyp)
                         DO i = 1,3
                            aveccof(i,ie,lm,na)=aveccof(i,ie,lm,na) +fgp(i)*ctmp*alo1(lo,ntyp)
                            bveccof(i,ie,lm,na)=bveccof(i,ie,lm,na) +fgp(i)*ctmp*blo1(lo,ntyp)
                            cveccof(i,m,ie,lo,na)=cveccof(i,m,ie,lo,na)+fgp(i)*ctmp*clo1(lo,ntyp)
                         ENDDO
160
                         IF (noco%l_soc.AND.sym%invs) THEN
161
                            ctmp = zMat%data_c(nbasf,ie) * CONJG(clotmp(m))*(-1)**(l-m)
162 163 164 165 166 167 168 169 170 171 172 173 174
                            na2 = sym%invsatnr(na)
                            lmp = ll1 - m
                            acof(ie,lmp,na2) = acof(ie,lmp,na2) +ctmp*alo1(lo,ntyp)
                            bcof(ie,lmp,na2) = bcof(ie,lmp,na2) +ctmp*blo1(lo,ntyp)
                            ccof(-m,ie,lo,na2) = ccof(-m,ie,lo,na2) + ctmp*clo1(lo,ntyp)
                            acoflo(-m,ie,lo,na2) = acoflo(-m,ie,lo,na2) +ctmp*alo1(lo,ntyp)
                            bcoflo(-m,ie,lo,na2) = bcoflo(-m,ie,lo,na2) +ctmp*blo1(lo,ntyp)
                            DO i = 1,3
                               aveccof(i,ie,lmp,na2)=aveccof(i,ie,lmp,na2)-fgp(i)*ctmp*alo1(lo,ntyp)
                               bveccof(i,ie,lmp,na2)=bveccof(i,ie,lmp,na2)-fgp(i)*ctmp*blo1(lo,ntyp)
                               cveccof(i,-m,ie,lo,na2) =cveccof(i,-m,ie,lo,na2) -fgp(i)*ctmp*clo1(lo,ntyp)
                            ENDDO
                         ENDIF
175 176 177 178 179 180 181
                      ENDDO ! loop over m
                   ENDDO    ! loop over eigenstates (ie)
                ELSE
                   nkvec(lo,na) = nkvec(lo,na) - 1
                END IF   ! linind
             END IF      ! nkvec(lo,na) < 2*(2*atoms%llo + 1)
          ELSE
182
             CALL juDFT_error("invsat =/= 0 or 1",calledby ="abclocdn_pulay")
183 184 185 186 187 188
          ENDIF
       ELSE
          enough(na) = .FALSE.
       ENDIF ! s > eps  & l >= 1  
    END DO
    IF ((k.EQ.nv) .AND. (.NOT.enough(na))) THEN
189 190
       WRITE (6,FMT=*) 'abclocdn did not find enough linearly independent'
       WRITE (6,FMT=*) 'ccof coefficient-vectors. the linear independence'
191 192
       WRITE (6,FMT=*) 'quality, linindq, is set to: ',linindq,'.'
       WRITE (6,FMT=*) 'this value might be to large.'
193
       CALL juDFT_error("did not find enough lin. ind. ccof-vectors" ,calledby ="abclocdn_pulay")
194 195 196 197
    END IF

  END SUBROUTINE abclocdn_pulay
END MODULE m_abclocdnpulay