File: genlmirep.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (86 lines) | stat: -rw-r--r-- 2,278 bytes parent folder | download
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
43
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

! Copyright (C) 2002-2008 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine genlmirep(lmax,ld,elm,ulm)
use modmain
implicit none
! arguments
integer, intent(in) :: lmax
integer, intent(in) :: ld
real(8), intent(out) :: elm(ld,natmtot)
complex(8), intent(out) :: ulm(ld,ld,natmtot)
! local variables
integer isym,lspl,is,ia,ias
integer lmmax,i,j,l,lm,n,p
integer info,lwork
! allocatable arrays
real(8), allocatable :: rwork(:)
complex(8), allocatable :: ulat(:,:,:)
complex(8), allocatable :: a(:,:),b(:,:)
complex(8), allocatable :: h(:,:)
complex(8), allocatable :: work(:)
lmmax=(lmax+1)**2
allocate(rwork(3*lmmax))
allocate(ulat(lmmax,lmmax,nsymlat))
allocate(a(lmmax,lmmax),b(lmmax,lmmax))
allocate(h(lmmax,lmmax))
lwork=2*lmmax
allocate(work(lwork))
! construct (l,m) rotation matrix for each lattice symmetry
a(:,:)=0.d0
do i=1,lmmax
  a(i,i)=1.d0
end do
do isym=1,nsymlat
  call rotzflm(symlatc(:,:,isym),0,lmax,lmmax,lmmax,lmmax,a,ulat(:,:,isym))
end do
! set up pseudorandom symmetric matrix H
h(:,:)=0.d0
p=1
do l=0,lmax
  n=2*l+1
  lm=idxlm(l,-l)
  do i=lm,lm+n-1
    do j=i,lm+n-1
! Park and Miller linear congruential generator
      p=mod(p*171,30269)
      h(i,j)=mod(p,lmmax)
      h(j,i)=h(i,j)
    end do
  end do
end do
! loop over species and atoms
do is=1,nspecies
  do ia=1,natoms(is)
    ias=idxas(ia,is)
! symmetrise H with site symmetries
    b(:,:)=0.d0
    do isym=1,nsymsite(ias)
! spatial rotation element in lattice point group
      lspl=lsplsyms(isym,ias)
! apply lattice symmetry as U*H*conjg(U')
      call zgemm('N','N',lmmax,lmmax,lmmax,zone,ulat(:,:,lspl),lmmax,h,lmmax, &
       zzero,a,lmmax)
      call zgemm('N','C',lmmax,lmmax,lmmax,zone,a,lmmax,ulat(:,:,lspl),lmmax, &
       zone,b,lmmax)
    end do
! block diagonalise symmetrised H
    do l=0,lmax
      n=2*l+1
      lm=idxlm(l,-l)
      call zheev('V','U',n,b(lm,lm),lmmax,elm(lm,ias),work,lwork,rwork,info)
    end do
! the unitary matrix U is the transpose of the eigenvector array
    do i=1,lmmax
      do j=1,lmmax
        ulm(i,j,ias)=b(j,i)
      end do
    end do
  end do
end do
deallocate(rwork,ulat,a,b,h,work)
return
end subroutine