hf_setup.F90 12.6 KB
Newer Older
1

Daniel Wortmann's avatar
Daniel Wortmann committed
2
MODULE m_hf_setup
Daniel Wortmann's avatar
Daniel Wortmann committed
3
CONTAINS
Daniel Wortmann's avatar
Daniel Wortmann committed
4
  SUBROUTINE hf_setup(hybrid,input,sym,kpts,DIMENSION,atoms,mpi,noco,cell,oneD,results,jsp,eig_id_hf,&
Daniel Wortmann's avatar
Daniel Wortmann committed
5
       hybdat,irank2,it,l_real,vr0,eig_irr) 
Daniel Wortmann's avatar
Daniel Wortmann committed
6
7
8
9
10
11
12
13
14
    USE m_types
    USE m_eig66_io
    USE m_util
    USE m_apws
    USE m_checkolap
    USE m_read_core
    USE m_gen_wavf
    IMPLICIT NONE
    TYPE(t_hybrid),INTENT(INOUT):: hybrid
15
    TYPE(t_kpts),INTENT(IN)     :: kpts
Daniel Wortmann's avatar
Daniel Wortmann committed
16
    TYPE(t_dimension),INTENT(IN):: DIMENSION
17
18
19
20
21
22
23
    TYPE(t_atoms),INTENT(IN)    :: atoms
    TYPE(t_mpi),INTENT(IN)      :: mpi
    TYPE(t_noco),INTENT(IN)     :: noco
    TYPE(t_cell),INTENT(IN)     :: cell
    TYPE(t_oneD),INTENT(IN)     :: oneD
    TYPE(t_input),INTENT(IN)    :: input
    TYPE(t_sym),INTENT(IN)      :: sym
Daniel Wortmann's avatar
Daniel Wortmann committed
24
    TYPE(t_results),INTENT(INOUT):: results
25
    INTEGER,INTENT(IN)          :: irank2(:),it
Daniel Wortmann's avatar
Daniel Wortmann committed
26

27
28
    INTEGER,INTENT(IN)          :: jsp,eig_id_hf
    REAL, INTENT(IN)            :: vr0(:,:,:)
Daniel Wortmann's avatar
Daniel Wortmann committed
29
30
    LOGICAL,INTENT(IN)          :: l_real
    TYPE(t_hybdat),INTENT(INOUT):: hybdat
Daniel Wortmann's avatar
Daniel Wortmann committed
31
32
    REAL,ALLOCATABLE,INTENT(OUT):: eig_irr(:,:)
 
33

Daniel Wortmann's avatar
Daniel Wortmann committed
34
    INTEGER:: ok,nk,nrec1,i,j,ll,l1,l2,ng,itype,n,l,n1,n2,nn
Daniel Wortmann's avatar
Daniel Wortmann committed
35
36
37
38
39
40
41
42


    TYPE(t_zmat),ALLOCATABLE :: zmat(:)
    REAL,    ALLOCATABLE    ::  basprod(:)
    REAL                    ::  el_eig(0:atoms%lmaxd,atoms%ntype), ello_eig(atoms%nlod,atoms%ntype),bk(3)
    INTEGER                 ::  degenerat(DIMENSION%neigd2+1,kpts%nkpt)
    INTEGER                 :: matind(DIMENSION%nbasfcn,2),nred
    TYPE(t_lapw)            :: lapw
Daniel Wortmann's avatar
Daniel Wortmann committed
43
  
Daniel Wortmann's avatar
Daniel Wortmann committed
44
45
46
47
    LOGICAL:: skip_kpt(kpts%nkpt)
    REAL   :: g(3)
    skip_kpt=.FALSE.

48
    IF( hybrid%l_calhf ) THEN
Daniel Wortmann's avatar
Daniel Wortmann committed
49
       !
50
       !             preparations for HF and hybrid functional calculation
Daniel Wortmann's avatar
Daniel Wortmann committed
51
       !
52
       CALL timestart("gen_bz and gen_wavf")
Daniel Wortmann's avatar
Daniel Wortmann committed
53

Daniel Wortmann's avatar
Daniel Wortmann committed
54
       ALLOCATE(zmat(kpts%nkptf),stat=ok)
Daniel Wortmann's avatar
Daniel Wortmann committed
55
56
       IF( ok .NE. 0 ) STOP 'eigen_hf: failure allocation z_c'
       ALLOCATE ( eig_irr(DIMENSION%neigd2,kpts%nkpt)      ,stat=ok )
Daniel Wortmann's avatar
Daniel Wortmann committed
57
       IF( ok .NE. 0 ) STOP 'eigen_hf: failure allocation eig_irr'
58
       ALLOCATE ( hybdat%kveclo_eig(atoms%nlotot,kpts%nkpt)  ,stat=ok )
Daniel Wortmann's avatar
Daniel Wortmann committed
59
       IF( ok .NE. 0 ) STOP 'eigen_hf: failure allocation hybdat%kveclo_eig'
60
61
       eig_irr = 0 ; hybdat%kveclo_eig = 0

Daniel Wortmann's avatar
Daniel Wortmann committed
62
63
64

       zmat(:)%l_real=l_real
      ! Reading the eig file
Daniel Wortmann's avatar
Daniel Wortmann committed
65
       DO nk = 1,kpts%nkpt
66
#               ifdef CPP_MPI
Daniel Wortmann's avatar
Daniel Wortmann committed
67
68
          ! jump to next k-point if this process is not present in communicator
          IF ( skip_kpt(nk) ) CYCLE
69
#               endif
Daniel Wortmann's avatar
Daniel Wortmann committed
70
71
72
73
74
          nrec1 = kpts%nkpt*(jsp-1) + nk
          zmat(nk)%nbasfcn=dimension%nbasfcn
          zmat(nk)%nbands=dimension%neigd2
          if (l_real) THEN
             ALLOCATE(zmat(nk)%z_r(dimension%nbasfcn,dimension%neigd2))
Daniel Wortmann's avatar
Daniel Wortmann committed
75
             ALLOCATE(zmat(nk)%z_c(0,0))
Daniel Wortmann's avatar
Daniel Wortmann committed
76
77
          else
             ALLOCATE(zmat(nk)%z_c(dimension%nbasfcn,dimension%neigd2))
Daniel Wortmann's avatar
Daniel Wortmann committed
78
             ALLOCATE(zmat(nk)%z_r(0,0))
Daniel Wortmann's avatar
Daniel Wortmann committed
79
          endif
Daniel Wortmann's avatar
Daniel Wortmann committed
80
          CALL read_eig(eig_id_hf,nk,jsp,el=el_eig,ello=ello_eig, neig=hybrid%ne_eig(nk),eig=eig_irr(:,nk), w_iks=results%w_iks(:,nk,jsp),kveclo=hybdat%kveclo_eig(:,nk),zmat=zmat(nk))
81
     
Daniel Wortmann's avatar
Daniel Wortmann committed
82
       END DO
Daniel Wortmann's avatar
Daniel Wortmann committed
83
84
85
86
87
88
89
90
91
92
       !Allocate further space
       DO nk=kpts%nkpt+1,kpts%nkptf
          zmat(nk)%nbasfcn=dimension%nbasfcn
          zmat(nk)%nbands=dimension%neigd2
          if (l_real) THEN
             ALLOCATE(zmat(nk)%z_r(dimension%nbasfcn,dimension%neigd2))
          else
             ALLOCATE(zmat(nk)%z_c(dimension%nbasfcn,dimension%neigd2))
          endif
       Enddo
Daniel Wortmann's avatar
Daniel Wortmann committed
93
94
95
96
97
98
99
100
101
102
103
104
105
       !
       !determine degenerate states at each k-point
       !
       ! degenerat(i) =1  band i  is not degenerat ,
       ! degenerat(i) =j  band i  has j-1 degenart states ( i, i+1, ..., i+j)
       ! degenerat(i) =0  band i  is  degenerat, but is not the lowest band
       !                  of the group of degenerate states
       IF ( mpi%irank == 0 ) THEN
          WRITE(6,*)
          WRITE(6,'(A)') "   k-point      |   number of occupied&
               &bands  |   maximal number of bands"
       END IF
       degenerat = 1
Daniel Wortmann's avatar
Daniel Wortmann committed
106
       hybrid%nobd      = 0
Daniel Wortmann's avatar
Daniel Wortmann committed
107
       DO nk=1 ,kpts%nkpt
108
#               ifdef CPP_MPI
Daniel Wortmann's avatar
Daniel Wortmann committed
109
110
          ! jump to next k-point if this k-point is not treated at this process
          IF ( skip_kpt(nk) ) CYCLE
111
#               endif
Daniel Wortmann's avatar
Daniel Wortmann committed
112
113
          DO i=1,hybrid%ne_eig(nk)
             DO j=i+1,hybrid%ne_eig(nk)
Daniel Wortmann's avatar
Daniel Wortmann committed
114
115
                IF( ABS(eig_irr(i,nk)-eig_irr(j,nk)) < 1E-07) THEN !0.015
                   degenerat(i,nk) = degenerat(i,nk) + 1
116
                END IF
Daniel Wortmann's avatar
Daniel Wortmann committed
117
118
119
             END DO
          END DO

Daniel Wortmann's avatar
Daniel Wortmann committed
120
          DO i=1,hybrid%ne_eig(nk)
Daniel Wortmann's avatar
Daniel Wortmann committed
121
122
123
124
125
126
127
             IF( degenerat(i,nk) .NE. 1 .OR. degenerat(i,nk) .NE. 0 ) &
                  degenerat(i+1:i+degenerat(i,nk)-1,nk) = 0
          END DO


          ! set the size of the exchange matrix in the space of the wavefunctions

Daniel Wortmann's avatar
Daniel Wortmann committed
128
129
          hybrid%nbands(nk)=hybrid%bands1
          IF(hybrid%nbands(nk).GT.hybrid%ne_eig(nk)) THEN
Daniel Wortmann's avatar
Daniel Wortmann committed
130
             IF ( mpi%irank == 0 ) THEN
Daniel Wortmann's avatar
Daniel Wortmann committed
131
                WRITE(*,*) ' maximum for hybrid%nbands is', hybrid%ne_eig(nk)
Daniel Wortmann's avatar
Daniel Wortmann committed
132
                WRITE(*,*) ' increase energy window to obtain enough eigenvalues'
Daniel Wortmann's avatar
Daniel Wortmann committed
133
                WRITE(*,*) ' set hybrid%nbands equal to hybrid%ne_eig'
Daniel Wortmann's avatar
Daniel Wortmann committed
134
             END IF
Daniel Wortmann's avatar
Daniel Wortmann committed
135
             hybrid%nbands(nk)=hybrid%ne_eig(nk)
Daniel Wortmann's avatar
Daniel Wortmann committed
136
          END IF
137

Daniel Wortmann's avatar
Daniel Wortmann committed
138
139
140
          DO i = hybrid%nbands(nk)-1,1,-1
             IF( (degenerat(i,nk) .GE. 1) .AND. (degenerat(i,nk)+i-1 .NE. hybrid%nbands(nk) ) ) THEN
                hybrid%nbands(nk) = i + degenerat(i,nk) - 1
Daniel Wortmann's avatar
Daniel Wortmann committed
141
142
143
                EXIT
             END IF
          END DO
144

Daniel Wortmann's avatar
Daniel Wortmann committed
145
146
          DO i = 1,hybrid%ne_eig(nk)
             IF(results%w_iks(i,nk,jsp) .GT. 0d0 ) hybrid%nobd(nk) = hybrid%nobd(nk) + 1
147

Daniel Wortmann's avatar
Daniel Wortmann committed
148
          END DO
Daniel Wortmann's avatar
Daniel Wortmann committed
149
          IF (hybrid%nobd(nk)>hybrid%nbands(nk)) THEN
Daniel Wortmann's avatar
Daniel Wortmann committed
150
             CALL judft_warn("More occupied bands than total no of bands!?")
Daniel Wortmann's avatar
Daniel Wortmann committed
151
             hybrid%nbands(nk)=hybrid%nobd(nk)
Daniel Wortmann's avatar
Daniel Wortmann committed
152
          ENDIF
Daniel Wortmann's avatar
Daniel Wortmann committed
153
          PRINT *,"bands:",nk, hybrid%nobd(nk),hybrid%nbands(nk),hybrid%ne_eig(nk)
Daniel Wortmann's avatar
Daniel Wortmann committed
154
       END DO
155
156

#             ifdef CPP_MPI
Daniel Wortmann's avatar
Daniel Wortmann committed
157
158
159
160
161
       ! send results for occupied bands to all processes
       sndreqd = 0 ; rcvreqd = 0
       DO nk = 1,kpts%nkpt
          IF ( skip_kpt(nk) ) THEN
             rcvreqd = rcvreqd + 1
Daniel Wortmann's avatar
Daniel Wortmann committed
162
             CALL MPI_IRECV(hybrid%nobd(nk),1,MPI_INTEGER4, MPI_ANY_SOURCE,TAG_SNDRCV_HYBDAT%NOBD+nk, mpi,rcvreq(rcvreqd),ierr(1))
Daniel Wortmann's avatar
Daniel Wortmann committed
163
164
165
166
          ELSE
             i = MOD( mpi%irank + isize2(nk), mpi%isize )
             DO WHILE ( i < mpi%irank-irank2(nk) .OR. i >= mpi%irank-irank2(nk)+isize2(nk) )
                sndreqd = sndreqd + 1
Daniel Wortmann's avatar
Daniel Wortmann committed
167
                CALL MPI_ISSEND(hybrid%nobd(nk),1,MPI_INTEGER4,i, TAG_SNDRCV_HYBDAT%NOBD+nk,mpi, sndreq(sndreqd),ierr(1) )
Daniel Wortmann's avatar
Daniel Wortmann committed
168
169
170
171
172
173
                i = MOD( i + isize2(nk), mpi%isize )
             END DO
          END IF
       END DO
       CALL MPI_WAITALL( rcvreqd, rcvreq, MPI_STATUSES_IGNORE, ierr(1) )
       ! Necessary to avoid compiler optimization
Daniel Wortmann's avatar
Daniel Wortmann committed
174
175
       ! Compiler does not know that hybrid%nobd is modified in mpi_waitall
       CALL MPI_GET_ADDRESS( hybrid%nobd, addr, ierr(1) )
Daniel Wortmann's avatar
Daniel Wortmann committed
176
       rcvreqd = 0
177
178
179

#             endif

Daniel Wortmann's avatar
Daniel Wortmann committed
180
       ! spread hybrid%nobd from IBZ to whole BZ
Daniel Wortmann's avatar
Daniel Wortmann committed
181
182
       DO nk = 1,kpts%nkptf
          i       = kpts%bkp(nk)
Daniel Wortmann's avatar
Daniel Wortmann committed
183
          hybrid%nobd(nk)= hybrid%nobd(i)
Daniel Wortmann's avatar
Daniel Wortmann committed
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
       END DO

       !
       ! generate eigenvectors z and MT coefficients from the previous iteration
       ! at all k-points
       !
       CALL gen_wavf(&
            &                 kpts%nkpt,kpts,it,sym,&
            &                 atoms,el_eig,&
            &                 ello_eig,cell,DIMENSION,&
            &                 hybrid,vr0,&
            &                 hybdat,&
            &                 noco,oneD,mpi,irank2,&
            &                 input,jsp,&
            &                 zmat)

       ! generate core wave functions (-> core1/2(jmtd,hybdat%nindxc,0:lmaxc,ntype) )
       CALL corewf(atoms,jsp,input,DIMENSION,vr0,&
            &                    hybdat%lmaxcd,hybdat%maxindxc,mpi,&
            &                    hybdat%lmaxc,hybdat%nindxc,hybdat%core1,hybdat%core2,hybdat%eig_c)
204
205

#             ifdef CPP_MPI
Daniel Wortmann's avatar
Daniel Wortmann committed
206
207
       ! wait until all files are written in gen_wavf
       CALL MPI_BARRIER(mpi%mpi_comm,ierr)
208
209
#             endif

Daniel Wortmann's avatar
Daniel Wortmann committed
210
211
212
213
214
215
216
217
218
219
220
221
222
223
       !
       ! check olap between core-basis/core-valence/basis-basis
       !
       CALL checkolap(atoms,hybdat,hybrid,&
            &                      kpts%nkpt,kpts,&
            &                      DIMENSION,mpi,irank2,skip_kpt,&
            &                      input,sym,&
            &                      noco,cell,lapw,jsp)

       !
       ! set up pointer pntgpt
       !

       ! setup dimension of pntgpt
Daniel Wortmann's avatar
Daniel Wortmann committed
224
       ALLOCATE(hybdat%pntgptd(3))
Daniel Wortmann's avatar
Daniel Wortmann committed
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
       hybdat%pntgptd = 0
       DO nk = 1,kpts%nkptf
          CALL apws(DIMENSION,input,noco, kpts,nk,cell,sym%zrfs,&
               &                    1,jsp, bk,lapw,matind,nred)
          hybdat%pntgptd(1) = MAXVAL( (/ ( ABS(lapw%k1(i,jsp)),i=1,lapw%nv(jsp)), hybdat%pntgptd(1) /) )
          hybdat%pntgptd(2) = MAXVAL( (/ ( ABS(lapw%k2(i,jsp)),i=1,lapw%nv(jsp)), hybdat%pntgptd(2) /) )
          hybdat%pntgptd(3) = MAXVAL( (/ ( ABS(lapw%k3(i,jsp)),i=1,lapw%nv(jsp)), hybdat%pntgptd(3) /) )
       END DO

       ALLOCATE( hybdat%pntgpt(-hybdat%pntgptd(1):hybdat%pntgptd(1), -hybdat%pntgptd(2):hybdat%pntgptd(2),&
            &                         -hybdat%pntgptd(3):hybdat%pntgptd(3),kpts%nkptf),stat=ok )
       IF( ok .NE. 0 ) STOP 'eigen_hf: failure allocation pntgpt'
       hybdat%pntgpt = 0
       DO nk = 1,kpts%nkptf
          CALL apws( DIMENSION,input,noco, kpts,nk,cell,sym%zrfs,&
               &                     1,jsp, bk,lapw,matind,nred)
          DO i = 1,lapw%nv(jsp)
             g = (/ lapw%k1(i,jsp),lapw%k2(i,jsp),lapw%k3(i,jsp) /)
             hybdat%pntgpt(g(1),g(2),g(3),nk) = i
          END DO
       END DO

       ALLOCATE ( basprod(atoms%jmtd),stat=ok )
       IF( ok .NE. 0 )STOP 'eigen_hf: failure allocation basprod'
       ALLOCATE ( hybdat%prodm(hybrid%maxindxm1,hybrid%maxindxp1,0:hybrid%maxlcutm1,atoms%ntype), stat= ok )
       IF( ok .NE. 0 ) STOP 'eigen_hf: failure allocation hybdat%prodm'
       ALLOCATE ( hybdat%prod(hybrid%maxindxp1,0:hybrid%maxlcutm1,atoms%ntype),stat= ok )
       IF( ok .NE. 0 ) STOP 'eigen_hf: failure allocation hybdat%prod'
       basprod = 0 ; hybdat%prodm = 0 ; hybdat%prod%l1 = 0 ; hybdat%prod%l2 = 0
       hybdat%prod%n1 = 0 ; hybdat%prod%n2 = 0
255
256
       ALLOCATE(hybdat%nindxp1(0:hybrid%maxlcutm1,atoms%ntype))
       hybdat%nindxp1 = 0
Daniel Wortmann's avatar
Daniel Wortmann committed
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
       DO itype = 1,atoms%ntype
          ng = atoms%jri(itype)
          DO l2 = 0,MIN(atoms%lmax(itype),hybrid%lcutwf(itype))
             ll = l2
             DO l1 = 0,ll
                IF(ABS(l1-l2).LE.hybrid%lcutm1(itype)) THEN
                   DO n2 = 1,hybrid%nindx(l2,itype)
                      nn = hybrid%nindx(l1,itype)
                      IF(l1.EQ.l2) nn = n2
                      DO n1 = 1,nn
                         ! Calculate all basis-function hybdat%products to obtain
                         ! the overlaps with the hybdat%product-basis functions (hybdat%prodm)
                         basprod(:ng) = ( hybdat%bas1(:ng,n1,l1,itype)*hybdat%bas1(:ng,n2,l2,itype) +hybdat%bas2(:ng,n1,l1,itype)*hybdat%bas2(:ng,n2,l2,itype)) / atoms%rmsh(:ng,itype)
                         DO l = ABS(l1-l2),MIN(hybrid%lcutm1(itype),l1+l2)
                            IF(MOD(l1+l2+l,2).EQ.0) THEN
272
273
                               hybdat%nindxp1(l,itype)    = hybdat%nindxp1(l,itype) + 1
                               n                  = hybdat%nindxp1(l,itype)
Daniel Wortmann's avatar
Daniel Wortmann committed
274
275
276
277
278
279
280
                               hybdat%prod(n,l,itype)%l1 = l1
                               hybdat%prod(n,l,itype)%l2 = l2
                               hybdat%prod(n,l,itype)%n1 = n1
                               hybdat%prod(n,l,itype)%n2 = n2
                               DO i = 1,hybrid%nindxm1(l,itype)
                                  hybdat%prodm(i,n,l,itype) = intgrf(basprod(:ng)*hybrid%basm1(:ng,i,l,itype), atoms%jri,atoms%jmtd,atoms%rmsh,atoms%dx,atoms%ntype,itype,hybdat%gridf)
                               END DO
281
                            END IF
Daniel Wortmann's avatar
Daniel Wortmann committed
282
                         END DO
283
284

                      END DO
Daniel Wortmann's avatar
Daniel Wortmann committed
285
286
287
288
289
290
291
292
293
294
295
296
297
298
                   END DO
                END IF
             END DO
          END DO
       END DO
       DEALLOCATE(basprod)
       CALL timestop("gen_bz and gen_wavf")


    ELSE IF (hybrid%l_hybrid ) THEN ! hybrid%l_calhf is false

       ! Reading the eig file
       !DO nk = n_start,kpts%nkpt,n_stride
       DO nk = 1,kpts%nkpt,1
Daniel Wortmann's avatar
Daniel Wortmann committed
299
300
          CALL read_eig(eig_id_hf,nk,jsp,el=el_eig, ello=ello_eig,neig=hybrid%ne_eig(nk),w_iks=results%w_iks(:,nk,jsp))
          hybrid%nobd(nk) = COUNT(results%w_iks(:hybrid%ne_eig(nk),nk,jsp) > 0.0 )
Daniel Wortmann's avatar
Daniel Wortmann committed
301
       END DO
Daniel Wortmann's avatar
Daniel Wortmann committed
302
    
Daniel Wortmann's avatar
Daniel Wortmann committed
303
       hybrid%maxlmindx = MAXVAL((/ ( SUM( (/ (hybrid%nindx(l,itype)*(2*l+1), l=0,atoms%lmax(itype)) /) ),itype=1,atoms%ntype) /) )
Daniel Wortmann's avatar
Daniel Wortmann committed
304
       hybrid%nbands    = MIN( hybrid%bands1, DIMENSION%neigd )
Daniel Wortmann's avatar
Daniel Wortmann committed
305
306

    ENDIF ! hybrid%l_calhf
Daniel Wortmann's avatar
Daniel Wortmann committed
307
308
  END SUBROUTINE hf_setup
END MODULE m_hf_setup