File: potxculr.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (180 lines) | stat: -rw-r--r-- 4,838 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
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180

! Copyright (C) 2018 T. Mueller, J. K. Dewhurst, S. Sharma and E. K. U. Gross.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine potxculr
use modmain
use modulr
use modomp
implicit none
! local variables
integer idm,is,ias,ir
integer np,npc,i
integer nthd,ithd
real(8) t1
! allocatable arrays
real(8), allocatable :: vxcrmt(:,:,:),vxcrir(:,:)
real(8), allocatable :: bxcrmt(:,:,:,:),bxcrir(:,:,:)
real(8), allocatable :: rhomt_(:,:),magmt_(:,:,:)
real(8), allocatable :: vxcmt_(:,:),bxcmt_(:,:,:)
real(8), allocatable :: rfmt(:,:)
complex(8), allocatable :: zfft(:,:)
allocate(vxcrmt(npcmtmax,natmtot,nqpt),vxcrir(ngtot,nqpt))
if (spinpol) then
  allocate(bxcrmt(npcmtmax,natmtot,ndmag,nqpt))
  allocate(bxcrir(ngtot,ndmag,nqpt))
end if
! generate the core density in spherical coordinates
allocate(rfmt(npmtmax,natmtot))
do ias=1,natmtot
  is=idxis(ias)
  i=1
  do ir=1,nrmti(is)
    t1=rhocr(ir,ias,1)*y00
    rfmt(i:i+lmmaxi-1,ias)=t1
    i=i+lmmaxi
  end do
  do ir=nrmti(is)+1,nrmt(is)
    t1=rhocr(ir,ias,1)*y00
    rfmt(i:i+lmmaxo-1,ias)=t1
    i=i+lmmaxo
  end do
end do
call omp_hold(nqpt,nthd)
!$OMP PARALLEL DEFAULT(SHARED) &
!$OMP PRIVATE(rhomt_,vxcmt_,magmt_,bxcmt_) &
!$OMP PRIVATE(ias,is,np,idm) &
!$OMP NUM_THREADS(nthd)
!$OMP DO
do ir=1,nqpt
  allocate(rhomt_(npmtmax,natmtot),vxcmt_(npmtmax,natmtot))
  if (spinpol) then
    allocate(magmt_(npmtmax,natmtot,ndmag),bxcmt_(npmtmax,natmtot,ndmag))
  end if
! convert the density from a coarse to a fine radial mesh
  do ias=1,natmtot
    is=idxis(ias)
    call dcopy(npcmt(is),rhormt(:,ias,ir),1,rhomt_(:,ias),1)
  end do
  call rfmtctof(rhomt_)
! add the core density
  do ias=1,natmtot
    is=idxis(ias)
    np=npmt(is)
    rhomt_(1:np,ias)=rhomt_(1:np,ias)+rfmt(1:np,ias)
  end do
! convert magnetisation from a coarse to a fine radial mesh
  do idm=1,ndmag
    do ias=1,natmtot
      is=idxis(ias)
      call dcopy(npcmt(is),magrmt(:,ias,idm,ir),1,magmt_(:,ias,idm),1)
    end do
    call rfmtctof(magmt_(:,:,idm))
  end do
! calculate the exchange-correlation potential and magnetic field
  call potxc(.false.,xctype,rhomt_,rhorir(:,ir),magmt_,magrir(:,:,ir),taumt, &
   tauir,exmt,exir,ecmt,ecir,vxcmt_,vxcrir(:,ir),bxcmt_,bxcrir(:,:,ir),wxcmt, &
   wxcir)
! convert muffin-tin potential and field from fine to coarse radial mesh
  do ias=1,natmtot
    is=idxis(ias)
    call rfmtftoc(nrmt(is),nrmti(is),vxcmt_(:,ias),vxcrmt(:,ias,ir))
  end do
  do idm=1,ndmag
    do ias=1,natmtot
      is=idxis(ias)
      call rfmtftoc(nrmt(is),nrmti(is),bxcmt_(:,ias,idm),bxcrmt(:,ias,idm,ir))
    end do
  end do
  deallocate(rhomt_,vxcmt_)
  if (spinpol) deallocate(magmt_,bxcmt_)
end do
!$OMP END DO
!$OMP END PARALLEL
call omp_free(nthd)
deallocate(rfmt)
! add the external Coulomb potential and Fourier transform to Q-space
do ias=1,natmtot
  is=idxis(ias)
  npc=npcmt(is)
  call omp_hold(npc,nthd)
  allocate(zfft(nqpt,0:nthd-1))
!$OMP PARALLEL DEFAULT(SHARED) &
!$OMP PRIVATE(ithd) &
!$OMP NUM_THREADS(nthd)
!$OMP DO
  do i=1,npc
    ithd=omp_get_thread_num()
    zfft(:,ithd)=vxcrmt(i,ias,:)+vclru(:)
    call zfftifc(3,ngridq,-1,zfft(:,ithd))
    vsqmt(i,ias,:)=vsqmt(i,ias,:)+zfft(:,ithd)
  end do
!$OMP END DO
!$OMP END PARALLEL
  deallocate(zfft)
  call omp_free(nthd)
end do
call omp_hold(ngtot,nthd)
allocate(zfft(nqpt,0:nthd-1))
!$OMP PARALLEL DEFAULT(SHARED) &
!$OMP PRIVATE(ithd) &
!$OMP NUM_THREADS(nthd)
!$OMP DO
do ir=1,ngtot
  ithd=omp_get_thread_num()
  zfft(:,ithd)=vxcrir(ir,:)+vclru(:)
  call zfftifc(3,ngridq,-1,zfft(:,ithd))
  vsqir(ir,:)=vsqir(ir,:)+zfft(:,ithd)
end do
!$OMP END DO
!$OMP END PARALLEL
deallocate(zfft)
call omp_free(nthd)
deallocate(vxcrmt,vxcrir)
if (.not.spinpol) return
! add the external magnetic fields and Fourier transform to Q-space
do idm=1,ndmag
  do ias=1,natmtot
    is=idxis(ias)
    npc=npcmt(is)
    call omp_hold(npc,nthd)
    allocate(zfft(nqpt,0:nthd-1))
!$OMP PARALLEL DEFAULT(SHARED) &
!$OMP PRIVATE(ithd) &
!$OMP NUM_THREADS(nthd)
!$OMP DO
    do i=1,npc
      ithd=omp_get_thread_num()
      zfft(:,ithd)=bxcrmt(i,ias,idm,:)+bfcru(:,idm)+bfcmtru(:,ias,idm)
      call zfftifc(3,ngridq,-1,zfft(:,ithd))
      bsqmt(i,ias,idm,:)=zfft(:,ithd)
    end do
!$OMP END DO
!$OMP END PARALLEL
    deallocate(zfft)
    call omp_free(nthd)
  end do
end do
do idm=1,ndmag
  call omp_hold(ngtot,nthd)
  allocate(zfft(nqpt,0:nthd-1))
!$OMP PARALLEL DEFAULT(SHARED) &
!$OMP PRIVATE(ithd) &
!$OMP NUM_THREADS(nthd)
!$OMP DO
  do ir=1,ngtot
    ithd=omp_get_thread_num()
    zfft(:,ithd)=bxcrir(ir,idm,:)+bfcru(:,idm)
    call zfftifc(3,ngridq,-1,zfft(:,ithd))
    bsqir(ir,idm,:)=zfft(:,ithd)
  end do
!$OMP END DO
!$OMP END PARALLEL
  deallocate(zfft)
  call omp_free(nthd)
end do
if (spinpol) deallocate(bxcrmt,bxcrir)
return
end subroutine