hsvac.F90 8.38 KB
Newer Older
Daniel Wortmann's avatar
Daniel Wortmann committed
1 2 3 4 5
!--------------------------------------------------------------------------------
! 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.
!--------------------------------------------------------------------------------
6
MODULE m_hsvac
7
  USE m_juDFT
8
CONTAINS
Daniel Wortmann's avatar
Daniel Wortmann committed
9 10 11 12
  !-----------------------------------------------------------
  !Calculate the vacuum contribution to the Hamiltonian and
  !Overlap matrix
  !-----------------------------------------------------------
13
  SUBROUTINE hsvac(&
Daniel Wortmann's avatar
Daniel Wortmann committed
14
       vacuum,stars,DIMENSION, atoms,mpi,jsp,input,v,evac,cell,&
15
       lapw,sym, noco,hmat,smat)
Daniel Wortmann's avatar
Daniel Wortmann committed
16
 
17 18 19 20

    USE m_vacfun
    USE m_types
    IMPLICIT NONE
21 22 23 24 25 26 27 28 29
    TYPE(t_dimension),INTENT(IN)  :: DIMENSION
    TYPE(t_input),INTENT(IN)      :: input
    TYPE(t_vacuum),INTENT(IN)     :: vacuum
    TYPE(t_noco),INTENT(IN)       :: noco
    TYPE(t_sym),INTENT(IN)        :: sym
    TYPE(t_stars),INTENT(IN)      :: stars
    TYPE(t_cell),INTENT(IN)       :: cell
    TYPE(t_atoms),INTENT(IN)      :: atoms
    TYPE(t_lapw),INTENT(IN)       :: lapw
Daniel Wortmann's avatar
Daniel Wortmann committed
30 31 32
    TYPE(t_mpi),INTENT(IN)        :: mpi
    TYPE(t_potden),INTENT(IN)     :: v
    CLASS(t_mat),INTENT(INOUT)     :: hmat(:,:),smat(:,:)
33 34
    !     ..
    !     .. Scalar Arguments ..
Daniel Wortmann's avatar
Daniel Wortmann committed
35
    INTEGER, INTENT (IN) :: jsp
36 37
    !     ..
    !     .. Array Arguments ..
38
    REAL,    INTENT (IN) :: evac(2,DIMENSION%jspd)
39 40 41
    !     ..
    !     .. Local Scalars ..
    COMPLEX hij,sij,apw_lo,c_1
Daniel Wortmann's avatar
Daniel Wortmann committed
42 43 44
    REAL d2,gz,sign,th,wronk
    INTEGER i,i2,ii,jj,ik,j,jk,k,jspin,ipot,ii0,i0
    INTEGER ivac,irec,imz,igvm2,igvm2i,s1,s2
45 46 47 48
    INTEGER jspin1,jspin2,jmax,jsp_start,jsp_end
    INTEGER i_start,nc,nc_0
    !     ..
    !     .. Local Arrays ..
Daniel Wortmann's avatar
Daniel Wortmann committed
49
    INTEGER:: nv2(DIMENSION%jspd)
50 51 52 53
    INTEGER kvac1(DIMENSION%nv2d,DIMENSION%jspd),kvac2(DIMENSION%nv2d,DIMENSION%jspd)
    INTEGER map2(DIMENSION%nvd,DIMENSION%jspd)
    COMPLEX tddv(DIMENSION%nv2d,DIMENSION%nv2d),tduv(DIMENSION%nv2d,DIMENSION%nv2d)
    COMPLEX tudv(DIMENSION%nv2d,DIMENSION%nv2d),tuuv(DIMENSION%nv2d,DIMENSION%nv2d)
54
    COMPLEX vxy_help(stars%ng2-1)
55 56 57 58
    COMPLEX a(DIMENSION%nvd,DIMENSION%jspd),b(DIMENSION%nvd,DIMENSION%jspd)
    REAL ddnv(DIMENSION%nv2d,DIMENSION%jspd),dudz(DIMENSION%nv2d,DIMENSION%jspd)
    REAL duz(DIMENSION%nv2d,DIMENSION%jspd), udz(DIMENSION%nv2d,DIMENSION%jspd)
    REAL uz(DIMENSION%nv2d,DIMENSION%jspd)
59 60
    !     ..

61
    d2 = SQRT(cell%omtil/cell%area)
62 63 64 65 66 67 68 69 70 71 72 73

    !--->    set up mapping function from 3d-->2d lapws

    DO jspin = 1,input%jspins
       nv2(jspin) = 0
       k_loop:DO  k = 1,lapw%nv(jspin)
          DO  j = 1,nv2(jspin)
             IF (lapw%k1(k,jspin).EQ.kvac1(j,jspin)&
                  .AND. lapw%k2(k,jspin).EQ.kvac2(j,jspin)) THEN
                map2(k,jspin) = j
                CYCLE k_loop
             END IF
74
          ENDDO
75
          nv2(jspin) = nv2(jspin) + 1
76
          IF (nv2(jspin)>DIMENSION%nv2d)  CALL juDFT_error("hsvac:dimension%nv2d",calledby ="hsvac")
77 78 79
          kvac1(nv2(jspin),jspin) = lapw%k1(k,jspin)
          kvac2(nv2(jspin),jspin) = lapw%k2(k,jspin)
          map2(k,jspin) = nv2(jspin)
80
       ENDDO k_loop
81 82 83
    ENDDO
    !--->    loop over the two vacuua (1: upper; 2: lower)
    DO ivac = 1,2
Daniel Wortmann's avatar
Daniel Wortmann committed
84 85 86 87 88 89 90 91 92
       sign = 3. - 2.*ivac !+/- 1
       DO jspin1=MERGE(1,jsp,noco%l_noco),MERGE(2,jsp,noco%l_noco) !loop over global spin
          s1=MIN(SIZE(hmat,1),jspin1) !in colinear case s1=1
          DO jspin2=MERGE(1,jsp,noco%l_noco),MERGE(2,jsp,noco%l_noco) !loop over global spin
             s2=MIN(SIZE(hmat,1),jspin2) !in colinear case s2=1
             ipot=3
             IF (jspin1==jspin2) ipot=jspin1
          !--->       get the wavefunctions and set up the tuuv, etc matrices          
             jspin=jsp
93
          CALL vacfun(&
94
               vacuum,DIMENSION,stars,&
95
               jsp,input,noco,ipot,&
Daniel Wortmann's avatar
Daniel Wortmann committed
96
               sym, cell,ivac,evac(1,1),lapw%bkpt,v%vacxy(:,:,ivac,ipot),v%vacz(:,:,:),kvac1,kvac2,nv2,&
97 98 99 100 101 102 103 104 105 106
               tuuv,tddv,tudv,tduv,uz,duz,udz,dudz,ddnv,wronk)
          !
          !--->       generate a and b coeffficients
          !
          IF (noco%l_noco) THEN
             DO jspin = 1,input%jspins
                DO k = 1,lapw%nv(jspin)
                   gz = sign*cell%bmat(3,3)*lapw%k3(k,jspin)
                   i2 = map2(k,jspin)
                   th = gz*cell%z1
Daniel Wortmann's avatar
Daniel Wortmann committed
107
                   c_1 = CMPLX( COS(th), SIN(th) )/ (d2*wronk)
108 109
                   a(k,jspin) = - c_1 * CMPLX(dudz(i2,jspin), gz*udz(i2,jspin) )
                   b(k,jspin) =   c_1 * CMPLX(duz(i2,jspin), gz* uz(i2,jspin) )
110 111 112 113
                ENDDO
             ENDDO
          ELSE
             DO k = 1,lapw%nv(jsp)
Daniel Wortmann's avatar
Daniel Wortmann committed
114
                gz = sign*cell%bmat(3,3)*lapw%gvec(3,k,jsp)
115 116
                i2 = map2(k,jsp)
                th = gz*cell%z1
Daniel Wortmann's avatar
Daniel Wortmann committed
117 118 119
                c_1 = CMPLX( COS(th), SIN(th) )/ (d2*wronk)
                a(k,jsp) = - c_1 * CMPLX(dudz(i2,jsp), gz*udz(i2,jsp) )
                b(k,jsp) =   c_1 * CMPLX(duz(i2,jsp), gz* uz(i2,jsp) )
120 121 122
             ENDDO
          ENDIF
          !--->       update hamiltonian and overlap matrices
Daniel Wortmann's avatar
Daniel Wortmann committed
123 124 125
          IF (jspin1==jspin2) THEN
             DO  i = mpi%n_rank+1,lapw%nv(jspin1),mpi%n_size
                i0=(i-1)/mpi%n_size+1 !local column index
126
                ik = map2(i,jspin1)
Daniel Wortmann's avatar
Daniel Wortmann committed
127
                DO j = 1,i - 1 !TODO check noco case
128
                   !--->             overlap: only  (g-g') parallel=0       '
129 130 131
                   IF (map2(j,jspin1).EQ.ik) THEN
                      sij = CONJG(a(i,jspin1))*a(j,jspin1) + &
                           CONJG(b(i,jspin1))*b(j,jspin1)*ddnv(ik,jspin1)
132
                      !+APW_LO
133
                      IF (input%l_useapw) THEN
134 135 136 137
                         apw_lo = CONJG(a(i,jspin1)*  uz(ik,jspin1) + b(i,jspin1)* udz(ik,jspin1) ) &
                              * (a(j,jspin1)* duz(ik,jspin1) + b(j,jspin1)*dudz(ik,jspin1) )&
                              +      (a(j,jspin1)*  uz(ik,jspin1) + b(j,jspin1)* udz(ik,jspin1) ) &
                              * CONJG(a(i,jspin1)* duz(ik,jspin1) + b(i,jspin1)*dudz(ik,jspin1) )
138
                         !            IF (i.lt.10) write (3,'(2i4,2f20.10)') i,j,apw_lo
Daniel Wortmann's avatar
Daniel Wortmann committed
139 140
                         IF (hmat(1,1)%l_real) THEN
                            hmat(s1,s2)%data_r(j,i0) = hmat(s1,s2)%data_r(j,i0) + 0.25 * REAL(apw_lo) 
141
                         ELSE 
Daniel Wortmann's avatar
Daniel Wortmann committed
142
                            hmat(s1,s2)%data_c(j,i0) = hmat(s1,s2)%data_c(j,i0) + 0.25 * apw_lo
143
                         ENDIF
144
                      ENDIF
Daniel Wortmann's avatar
Daniel Wortmann committed
145 146 147
                      !Overlapp Matrix
                      IF (hmat(1,1)%l_real) THEN
                         smat(s1,s2)%data_r(j,i0) = smat(s1,s2)%data_r(j,i0) + REAL(sij)
148
                      ELSE 
Daniel Wortmann's avatar
Daniel Wortmann committed
149
                         smat(s1,s2)%data_c(j,i0) = smat(s1,s2)%data_c(j,i0) + sij
150
                      ENDIF
151 152
                   END IF
                ENDDO
Daniel Wortmann's avatar
Daniel Wortmann committed
153
                !Diagonal term of Overlapp matrix, Hamiltonian later
154
                sij = CONJG(a(i,jspin1))*a(i,jspin1) + CONJG(b(i,jspin1))*b(i,jspin1)*ddnv(ik,jspin1)
Daniel Wortmann's avatar
Daniel Wortmann committed
155 156
                IF (hmat(1,1)%l_real) THEN
                   smat(s1,s2)%data_r(j,i0) = smat(s1,s2)%data_r(j,i0) + REAL(sij)
157
                ELSE
Daniel Wortmann's avatar
Daniel Wortmann committed
158
                   smat(s1,s2)%data_c(j,i0) = smat(s1,s2)%data_c(j,i0) + sij
159
                ENDIF
160
             ENDDO
161 162 163
          ENDIF

          !--->    hamiltonian update
Daniel Wortmann's avatar
Daniel Wortmann committed
164 165
          DO  i = mpi%n_rank+1,lapw%nv(jspin2),mpi%n_size
             i0=(i-1)/mpi%n_size+1 !local column index
166
             ik = map2(i,jspin1)
Daniel Wortmann's avatar
Daniel Wortmann committed
167
             DO j = 1,i
168
                jk = map2(j,jspin2)
169 170 171 172 173 174 175
                IF (jspin2>jspin1) THEN
                   hij = CONJG(CONJG(a(j,jspin2))* (tuuv(jk,ik)*a(i,jspin1) +tudv(jk,ik)*b(i,jspin1))&
                        + CONJG(b(j,jspin2))* (tddv(jk,ik)*b(i,jspin1) +tduv(jk,ik)*a(i,jspin1)))
                ELSE
                   hij = CONJG(a(i,jspin1))* (tuuv(ik,jk)*a(j,jspin2) +tudv(ik,jk)*b(j,jspin2))&
                        + CONJG(b(i,jspin1))* (tddv(ik,jk)*b(j,jspin2) +tduv(ik,jk)*a(j,jspin2))
                ENDIF
Daniel Wortmann's avatar
Daniel Wortmann committed
176 177
                IF (hmat(1,1)%l_real) THEN
                   hmat(s1,s2)%data_r(j,i0) = hmat(s1,s2)%data_r(j,i0) + REAL(hij)
178
                ELSE
179
                   hmat(s1,s2)%data_c(j,i0) = hmat(s1,s2)%data_c(j,i0) + hij
180
                ENDIF
181 182 183 184 185 186 187 188
             ENDDO
          ENDDO

          !--->    end of loop over different parts of the potential matrix
       ENDDO

       !---> end of loop over vacua
    ENDDO
Daniel Wortmann's avatar
Daniel Wortmann committed
189
 ENDDO
190 191 192

  END SUBROUTINE hsvac
END MODULE m_hsvac