gen_wavf.F90 12.3 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
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!     This module generates the cmt coefficients and eigenvectors z   !
!     at all kpoints nkpt from the irreducible kpoints nkpti          !
!     and writes them out in cmt and z, respectively.                 !
!                                                 M.Betzinger(09/07)  !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Matthias Redies's avatar
Matthias Redies committed
13

14
MODULE m_gen_wavf
15

16
CONTAINS
17

Matthias Redies's avatar
Matthias Redies committed
18 19
   SUBROUTINE gen_wavf(nkpti, kpts, it, sym, atoms, el_eig, ello_eig, cell, dimension, hybrid, vr0, &
                       hybdat, noco, oneD, mpi, input, jsp, zmat)
20 21

      ! nkpti      ::     number of irreducible k-points
Matthias Redies's avatar
Matthias Redies committed
22
      ! nkpt       ::     number of all k-points
23 24 25 26

      USE m_radfun
      USE m_radflo
      USE m_abcof
Matthias Redies's avatar
Matthias Redies committed
27 28
      USE m_trafo, ONLY: waveftrafo_genwavf
      USE m_util, ONLY: modulo1
29 30
      USE m_olap
      USE m_types
31
      USE m_hyb_abcrot
Daniel Wortmann's avatar
Daniel Wortmann committed
32
      USE m_io_hybrid
33

34
      IMPLICIT NONE
35

Matthias Redies's avatar
Matthias Redies committed
36 37
      TYPE(t_hybdat), INTENT(INOUT) :: hybdat
      TYPE(t_mpi), INTENT(IN)    :: mpi
38
      TYPE(t_dimension), INTENT(IN)    :: dimension
Matthias Redies's avatar
Matthias Redies committed
39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54
      TYPE(t_oneD), INTENT(IN)    :: oneD
      TYPE(t_hybrid), INTENT(IN)    :: hybrid
      TYPE(t_input), INTENT(IN)    :: input
      TYPE(t_noco), INTENT(IN)    :: noco
      TYPE(t_sym), INTENT(IN)    :: sym
      TYPE(t_cell), INTENT(IN)    :: cell
      TYPE(t_kpts), INTENT(IN)    :: kpts
      TYPE(t_atoms), INTENT(IN)    :: atoms
      TYPE(t_mat), INTENT(IN)    :: zmat(:) !for all kpoints

      INTEGER, INTENT(IN)    :: nkpti, it
      INTEGER, INTENT(IN)    :: jsp

      REAL, INTENT(IN)    :: vr0(:, :, :)!(jmtd,ntype,jspd)
      REAL, INTENT(IN)    :: el_eig(0:atoms%lmaxd, atoms%ntype)
      REAL, INTENT(IN)    :: ello_eig(atoms%nlod, atoms%ntype)
55 56

      ! local scalars
Matthias Redies's avatar
Matthias Redies committed
57
      INTEGER                 :: ilo, idum, m, irecl_cmt, irecl_z
58
      COMPLEX                 :: cdum
Daniel Wortmann's avatar
Daniel Wortmann committed
59
      TYPE(t_mat)             :: zhlp
Matthias Redies's avatar
Matthias Redies committed
60 61 62
      INTEGER                 :: ikpt0, ikpt, itype, iop, ispin, ieq, indx, iatom
      INTEGER                 :: i, j, l, ll, lm, ng, ok
      COMPLEX                 :: img = (0.0, 1.0)
63

Matthias Redies's avatar
Matthias Redies committed
64
      INTEGER                 :: nodem, noded
65
      REAL                    :: wronk
66

67 68
      INTEGER                 :: lower, upper
      LOGICAL                 :: found
69

70
      ! local arrays
Matthias Redies's avatar
Matthias Redies committed
71
      INTEGER                 :: rrot(3, 3, sym%nsym)
72
      INTEGER                 :: map_lo(atoms%nlod)
Matthias Redies's avatar
Matthias Redies committed
73 74 75 76 77 78 79 80 81 82 83 84
      INTEGER                 :: iarr(0:atoms%lmaxd, atoms%ntype)
      COMPLEX, ALLOCATABLE     :: acof(:, :, :), bcof(:, :, :), ccof(:, :, :, :)

      COMPLEX, ALLOCATABLE     :: cmt(:, :, :), cmthlp(:, :, :)

      REAL                    :: vr(atoms%jmtd, atoms%ntype, input%jspins)
      REAL, ALLOCATABLE        :: f(:, :, :), df(:, :, :)

      REAL                    :: flo(atoms%jmtd, 2, atoms%nlod)
      REAL                    :: uuilon(atoms%nlod, atoms%ntype), duilon(atoms%nlod, atoms%ntype)
      REAL                    :: ulouilopn(atoms%nlod, atoms%nlod, atoms%ntype)

85
      REAL                    :: bkpt(3)
86

87 88
!     local arrays for abcof1
!      COMPLEX                 ::  a(nvd,0:lmd,natd,nkpti),b(nvd,0:lmd,natd,nkpti)
89

90 91
      TYPE(t_lapw)  :: lapw(kpts%nkptf)
      TYPE(t_usdus) :: usdus
92

Matthias Redies's avatar
Matthias Redies committed
93 94
      CALL usdus%init(atoms, input%jspins)
      CALL zhlp%alloc(zmat(1)%l_real, zmat(1)%matsize1, zmat(1)%matsize2)
Daniel Wortmann's avatar
Daniel Wortmann committed
95

96
      ! setup rotations in reciprocal space
97
      DO iop = 1, sym%nsym
Matthias Redies's avatar
Matthias Redies committed
98 99
         IF (iop <= sym%nop) THEN
            rrot(:, :, iop) = transpose(sym%mrot(:, :, sym%invtab(iop)))
100
         ELSE
Matthias Redies's avatar
Matthias Redies committed
101
            rrot(:, :, iop) = -rrot(:, :, iop - sym%nop)
102
         END IF
103 104 105 106
      END DO

      ! generate G-vectors, which fulfill |k+G|<rkmax
      ! for all k-points
107
      DO ikpt = 1, kpts%nkptf
Matthias Redies's avatar
Matthias Redies committed
108
         CALL lapw(ikpt)%init(input, noco, kpts, atoms, sym, ikpt, cell, sym%zrfs)
109 110
      END DO

111
      ! set spherical component of the potential from the previous iteration vr
112
      vr = vr0
113 114 115 116 117 118 119 120 121 122 123

!       ALLOCATE ( z_out(nbasfcn,neigd,nkpti),stat=ok )
!       IF ( ok .ne. 0) STOP 'gen_wavf: failure allocation z'
!       z_out = 0
!       z_out(:,:,:nkpti) = z_in

      ! calculate radial basis functions belonging to the
      ! potential vr stored in bas1 and bas2
      ! bas1 denotes the large component
      ! bas2    "     "  small component

Matthias Redies's avatar
Matthias Redies committed
124
      ALLOCATE (f(atoms%jmtd, 2, 0:atoms%lmaxd), df(atoms%jmtd, 2, 0:atoms%lmaxd))
125 126
      f = 0
      df = 0
127
      iarr = 2
128
      DO itype = 1, atoms%ntype
Matthias Redies's avatar
Matthias Redies committed
129
         IF (mpi%irank == 0) WRITE (6, FMT=8000) itype
130
         ng = atoms%jri(itype)
Matthias Redies's avatar
Matthias Redies committed
131 132 133 134 135 136 137 138 139 140 141 142 143 144
         DO l = 0, atoms%lmax(itype)
            CALL radfun(l, itype, jsp, el_eig(l, itype), vr(:, itype, jsp), atoms, f(:, :, l), df(:, :, l), usdus, nodem, noded, wronk)
            IF (mpi%irank == 0) WRITE (6, FMT=8010) l, el_eig(l, itype), usdus%us(l, itype, jsp), usdus%dus(l, itype, jsp), nodem, &
               usdus%uds(l, itype, jsp), usdus%duds(l, itype, jsp), noded, usdus%ddn(l, itype, jsp), wronk

            hybdat%bas1(1:ng, 1, l, itype) = f(1:ng, 1, l)
            hybdat%bas2(1:ng, 1, l, itype) = f(1:ng, 2, l)
            hybdat%bas1(1:ng, 2, l, itype) = df(1:ng, 1, l)
            hybdat%bas2(1:ng, 2, l, itype) = df(1:ng, 2, l)

            hybdat%bas1_MT(1, l, itype) = usdus%us(l, itype, jsp)
            hybdat%drbas1_MT(1, l, itype) = usdus%dus(l, itype, jsp)
            hybdat%bas1_MT(2, l, itype) = usdus%uds(l, itype, jsp)
            hybdat%drbas1_MT(2, l, itype) = usdus%duds(l, itype, jsp)
145 146
         END DO

Matthias Redies's avatar
Matthias Redies committed
147 148
         IF (atoms%nlo(itype) >= 1) THEN
            CALL radflo(atoms, itype, jsp, ello_eig, vr(:, itype, jsp), f, df, mpi, usdus, uuilon, duilon, ulouilopn, flo)
149

Matthias Redies's avatar
Matthias Redies committed
150 151 152 153 154 155
            DO ilo = 1, atoms%nlo(itype)
               iarr(atoms%llo(ilo, itype), itype) = iarr(atoms%llo(ilo, itype), itype) + 1
               hybdat%bas1(1:ng, iarr(atoms%llo(ilo, itype), itype), atoms%llo(ilo, itype), itype) = flo(1:ng, 1, ilo)
               hybdat%bas2(1:ng, iarr(atoms%llo(ilo, itype), itype), atoms%llo(ilo, itype), itype) = flo(1:ng, 2, ilo)
               hybdat%bas1_MT(iarr(atoms%llo(ilo, itype), itype), atoms%llo(ilo, itype), itype) = usdus%ulos(ilo, itype, jsp)
               hybdat%drbas1_MT(iarr(atoms%llo(ilo, itype), itype), atoms%llo(ilo, itype), itype) = usdus%dulos(ilo, itype, jsp)
156 157
            END DO
         END IF
158
      END DO
Matthias Redies's avatar
Matthias Redies committed
159
      DEALLOCATE (f, df)
160

161 162
#if CPP_DEBUG
      ! consistency check
Matthias Redies's avatar
Matthias Redies committed
163
      IF (.not. all(iarr == hybrid%nindx)) STOP 'gen_wavf: counting error'
164 165
#endif

Matthias Redies's avatar
Matthias Redies committed
166 167 168 169
8000  FORMAT(1x, /, /, ' wavefunction parameters for atom type', i3, ':', /, t32, 'radial function', t79, &
             'energy derivative', /, t3, 'l', t8, 'energy', t26, 'value', t39, 'derivative', t53, &
             'nodes', t68, 'value', t81, 'derivative', t95, 'nodes', t107, 'norm', t119, 'wronskian')
8010  FORMAT(i3, f10.5, 2(5x, 1p, 2e16.7, i5), 1p, 2e16.7)
170 171 172 173

      ! determine boundaries for parallel calculations
      lower = 1
      upper = nkpti
174
      found = .true.
Matthias Redies's avatar
Matthias Redies committed
175
      
Matthias Redies's avatar
Matthias Redies committed
176
      IF (.NOT. found) THEN
177
         upper = 0
178 179 180 181 182 183
      END IF

      ! calculate wavefunction expansion in the the MT region
      ! (acof,bcof,ccof) and APW-basis coefficients
      ! (a,b,bascofold_lo) at irred. kpoints

Matthias Redies's avatar
Matthias Redies committed
184 185 186 187 188 189
      ALLOCATE (acof(dimension%neigd, 0:dimension%lmd, atoms%nat), stat=ok)
      IF (ok /= 0) STOP 'gen_wavf: failure allocation acof'
      ALLOCATE (bcof(dimension%neigd, 0:dimension%lmd, atoms%nat), stat=ok)
      IF (ok /= 0) STOP 'gen_wavf: failure allocation bcof'
      ALLOCATE (ccof(-atoms%llod:atoms%llod, dimension%neigd, atoms%nlod, atoms%nat), stat=ok)
      IF (ok /= 0) STOP 'gen_wavf: failure allocation ccof'
190
      ALLOCATE (cmt(dimension%neigd2, hybrid%maxlmindx, atoms%nat), stat=ok)
Matthias Redies's avatar
Matthias Redies committed
191
      IF (ok /= 0) STOP 'gen_wavf: Failure allocation cmt'
192
      ALLOCATE (cmthlp(dimension%neigd2, hybrid%maxlmindx, atoms%nat), stat=ok)
Matthias Redies's avatar
Matthias Redies committed
193
      IF (ok /= 0) STOP 'gen_wavf: failure allocation cmthlp'
194 195 196

      DO ikpt0 = lower, upper

197
         acof = 0; bcof = 0; ccof = 0
198

199 200 201
         ! abcof calculates the wavefunction coefficients
         ! stored in acof,bcof,ccof
         lapw(ikpt0)%nmat = lapw(ikpt0)%nv(jsp) + atoms%nlotot
Matthias Redies's avatar
Matthias Redies committed
202 203 204
         CALL abcof(input, atoms, sym, cell, lapw(ikpt0), hybrid%nbands(ikpt0), usdus, noco, jsp, &!hybdat%kveclo_eig(:,ikpt0),&
                    oneD, acof(:hybrid%nbands(ikpt0), :, :), bcof(:hybrid%nbands(ikpt0), :, :), &
                    ccof(:, :hybrid%nbands(ikpt0), :, :), zmat(ikpt0))
205 206

! call was ...
Matthias Redies's avatar
Matthias Redies committed
207 208 209 210 211 212 213 214
         ! gpt(1,:,:,ikpt0),gpt(2,:,:,ikpt0),&
         ! gpt(3,:,:,ikpt0),ngpt(:,ikpt0),&!k1hlp,k2hlp,k3hlp,nvhlp,&
         !    ngpt(jsp,ikpt0)+nbands(ikpt0),z(:,:,ikpt0),&!nvhlp(jsp)+ &
         !   &usdus,&
         !    noco,&
         !    jsp,kveclo_eig(:ikpt0),oneD,oneD,&
         !    acof(:nbands(ikpt0),:,:),&
         !    bcof(:nbands(ikpt0),:,:),ccof(:,:nbands(ikpt0),:,:) )
215

216 217
         ! MT wavefunction coefficients are calculated in a local coordinate system rotate them in the global one

Matthias Redies's avatar
Matthias Redies committed
218 219
         CALL hyb_abcrot(hybrid, atoms, hybrid%nbands(ikpt0), sym, cell, oneD, acof(:hybrid%nbands(ikpt0), :, :), &
                         bcof(:hybrid%nbands(ikpt0), :, :), ccof(:, :hybrid%nbands(ikpt0), :, :))
220 221 222 223 224 225 226 227 228 229 230 231

         ! decorate acof, bcof, ccof with coefficient i**l and store them
         ! in the field cmt(neigd,nkpt,maxlmindx,nat), i.e.
         ! where maxlmindx subsumes l,m and nindx

         cmt = 0
         iatom = 0
         DO itype = 1, atoms%ntype
            DO ieq = 1, atoms%neq(itype)
               iatom = iatom + 1
               indx = 0
               DO l = 0, atoms%lmax(itype)
Matthias Redies's avatar
Matthias Redies committed
232
                  ll = l*(l + 1)
233 234 235 236 237
                  cdum = img**l

                  ! determine number of local orbitals with quantum number l
                  ! map returns the number of the local orbital of quantum
                  ! number l in the list of all local orbitals of the atom type
Matthias Redies's avatar
Matthias Redies committed
238
                  idum = 0
239
                  map_lo = 0
Matthias Redies's avatar
Matthias Redies committed
240
                  IF (hybrid%nindx(l, itype) > 2) THEN
241
                     DO j = 1, atoms%nlo(itype)
Matthias Redies's avatar
Matthias Redies committed
242
                        IF (atoms%llo(j, itype) == l) THEN
243 244 245 246
                           idum = idum + 1
                           map_lo(idum) = j
                        END IF
                     END DO
247 248
                  END IF

249 250
                  DO M = -l, l
                     lm = ll + M
Matthias Redies's avatar
Matthias Redies committed
251
                     DO i = 1, hybrid%nindx(l, itype)
252
                        indx = indx + 1
Matthias Redies's avatar
Matthias Redies committed
253 254 255 256
                        IF (i == 1) THEN
                           cmt(:, indx, iatom) = cdum*acof(:, lm, iatom)
                        ELSE IF (i == 2) THEN
                           cmt(:, indx, iatom) = cdum*bcof(:, lm, iatom)
257 258
                        ELSE
                           idum = i - 2
Matthias Redies's avatar
Matthias Redies committed
259
                           cmt(:, indx, iatom) = cdum*ccof(M, :, map_lo(idum), iatom)
260 261 262 263 264 265
                        END IF
                     END DO
                  END DO
               END DO
            END DO
         END DO
266

267
         ! write cmt at irreducible k-points in direct-access file cmt
Matthias Redies's avatar
Matthias Redies committed
268 269 270
         CALL write_cmt(cmt, ikpt0)
         CALL zhlp%alloc(zmat(ikpt0)%l_real, zmat(ikpt0)%matsize1, zmat(ikpt0)%matsize2)

271 272 273 274 275
         IF (zhlp%l_real) THEN
            zhlp%data_r = zmat(ikpt0)%data_r
         ELSE
            zhlp%data_c = zmat(ikpt0)%data_c
         END IF
276
         CALL write_z(zhlp, kpts%nkptf*(jsp - 1) + ikpt0)
277 278 279

         ! generate wavefunctions coefficients at all k-points from
         ! irreducible k-points
280

281
         DO ikpt = 1, kpts%nkptf
Matthias Redies's avatar
Matthias Redies committed
282
            IF ((kpts%bkp(ikpt) == ikpt0) .AND. (ikpt0 /= ikpt)) THEN
283
               iop = kpts%bksym(ikpt)
Matthias Redies's avatar
Matthias Redies committed
284
               CALL waveftrafo_genwavf(cmthlp, zhlp%data_r, zhlp%data_c, cmt(:, :, :), zmat(1)%l_real, zmat(ikpt0)%data_r(:, :), &
285
                                       zmat(ikpt0)%data_c(:, :), ikpt0, iop, atoms, hybrid, kpts, sym, jsp, zmat(ikpt0)%matsize1, dimension, &
Matthias Redies's avatar
Matthias Redies committed
286
                                       hybrid%nbands(ikpt0), cell, lapw(ikpt0), lapw(ikpt), .true.)
287

Matthias Redies's avatar
Matthias Redies committed
288
               CALL write_cmt(cmthlp, ikpt)
289
               CALL write_z(zhlp, kpts%nkptf*(jsp - 1) + ikpt)
290 291
            END IF
         END DO  !ikpt
292 293
      END DO !ikpt0

Matthias Redies's avatar
Matthias Redies committed
294 295
      DEALLOCATE (acof, bcof, ccof)
      DEALLOCATE (cmt, cmthlp)
296

297
   END SUBROUTINE gen_wavf
298

299
END MODULE m_gen_wavf