File: ggair_sp_2b.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 (141 lines) | stat: -rw-r--r-- 3,882 bytes parent folder | download | duplicates (3)
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

! Copyright (C) 2009 T. McQueen and J. K. Dewhurst.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

!BOP
! !ROUTINE: ggair_sp_2b
! !INTERFACE:
subroutine ggair_sp_2b(g2up,g2dn,gvup,gvdn,vxup,vxdn,vcup,vcdn,dxdgu2,dxdgd2, &
 dxdgud,dcdgu2,dcdgd2,dcdgud)
! !USES:
use modmain
! !DESCRIPTION:
!   Post processing step of interstitial gradients for GGA type 2. See routine
!   {\tt ggamt\_sp\_2a} for full details.
!
! !REVISION HISTORY:
!   Created November 2009 (JKD and TMcQ)
!EOP
!BOC
implicit none
real(8), intent(in) :: g2up(ngtot)
real(8), intent(in) :: g2dn(ngtot)
real(8), intent(in) :: gvup(ngtot,3)
real(8), intent(in) :: gvdn(ngtot,3)
real(8), intent(inout) :: vxup(ngtot)
real(8), intent(inout) :: vxdn(ngtot)
real(8), intent(inout) :: vcup(ngtot)
real(8), intent(inout) :: vcdn(ngtot)
real(8), intent(in) :: dxdgu2(ngtot)
real(8), intent(in) :: dxdgd2(ngtot)
real(8), intent(in) :: dxdgud(ngtot)
real(8), intent(in) :: dcdgu2(ngtot)
real(8), intent(in) :: dcdgd2(ngtot)
real(8), intent(in) :: dcdgud(ngtot)
! local variables
integer ig,ifg,i
! allocatable arrays
real(8), allocatable :: rfir(:)
complex(8), allocatable :: zfft1(:),zfft2(:)
allocate(rfir(ngtot))
allocate(zfft1(ngtot),zfft2(ngtot))
!------------------!
!     exchange     !
!------------------!
! compute grad dxdgu2
zfft1(:)=dxdgu2(:)
call zfftifc(3,ngridg,-1,zfft1)
! (grad dxdgu2).(grad rhoup)
rfir(:)=0.d0
do i=1,3
  zfft2(:)=0.d0
  do ig=1,ngvec
    ifg=igfft(ig)
    zfft2(ifg)=vgc(i,ig)*cmplx(-aimag(zfft1(ifg)),dble(zfft1(ifg)),8)
  end do
  call zfftifc(3,ngridg,1,zfft2)
  rfir(:)=rfir(:)+dble(zfft2(:))*gvup(:,i)
end do
vxup(:)=vxup(:)-2.d0*(rfir(:)+dxdgu2(:)*g2up(:))-dxdgud(:)*g2dn(:)
! compute grad dxdgd2
zfft1(:)=dxdgd2(:)
call zfftifc(3,ngridg,-1,zfft1)
! (grad dxdgd2).(grad rhodn)
rfir(:)=0.d0
do i=1,3
  zfft2(:)=0.d0
  do ig=1,ngvec
    ifg=igfft(ig)
    zfft2(ifg)=vgc(i,ig)*cmplx(-aimag(zfft1(ifg)),dble(zfft1(ifg)),8)
  end do
  call zfftifc(3,ngridg,1,zfft2)
  rfir(:)=rfir(:)+dble(zfft2(:))*gvdn(:,i)
end do
vxdn(:)=vxdn(:)-2.d0*(rfir(:)+dxdgd2(:)*g2dn(:))-dxdgud(:)*g2up(:)
! compute grad dxdgud
zfft1(:)=dxdgud(:)
call zfftifc(3,ngridg,-1,zfft1)
! (grad dxdgud).(grad rhodn) and (grad dxdgud).(grad rhoup)
do i=1,3
  zfft2(:)=0.d0
  do ig=1,ngvec
    ifg=igfft(ig)
    zfft2(ifg)=vgc(i,ig)*cmplx(-aimag(zfft1(ifg)),dble(zfft1(ifg)),8)
  end do
  call zfftifc(3,ngridg,1,zfft2)
  vxup(:)=vxup(:)-dble(zfft2(:))*gvdn(:,i)
  vxdn(:)=vxdn(:)-dble(zfft2(:))*gvup(:,i)
end do
!---------------------!
!     correlation     !
!---------------------!
! compute grad dcdgu2
zfft1(:)=dcdgu2(:)
call zfftifc(3,ngridg,-1,zfft1)
! (grad dcdgu2).(grad rhoup)
rfir(:)=0.d0
do i=1,3
  zfft2(:)=0.d0
  do ig=1,ngvec
    ifg=igfft(ig)
    zfft2(ifg)=vgc(i,ig)*cmplx(-aimag(zfft1(ifg)),dble(zfft1(ifg)),8)
  end do
  call zfftifc(3,ngridg,1,zfft2)
  rfir(:)=rfir(:)+dble(zfft2(:))*gvup(:,i)
end do
vcup(:)=vcup(:)-2.d0*(rfir(:)+dcdgu2(:)*g2up(:))-dcdgud(:)*g2dn(:)
! compute grad dcdgd2
zfft1(:)=dcdgd2(:)
call zfftifc(3,ngridg,-1,zfft1)
! (grad dcdgd2).(grad rhodn)
rfir(:)=0.d0
do i=1,3
  zfft2(:)=0.d0
  do ig=1,ngvec
    ifg=igfft(ig)
    zfft2(ifg)=vgc(i,ig)*cmplx(-aimag(zfft1(ifg)),dble(zfft1(ifg)),8)
  end do
  call zfftifc(3,ngridg,1,zfft2)
  rfir(:)=rfir(:)+dble(zfft2(:))*gvdn(:,i)
end do
vcdn(:)=vcdn(:)-2.d0*(rfir(:)+dcdgd2(:)*g2dn(:))-dcdgud(:)*g2up(:)
! compute grad dcdgud
zfft1(:)=dcdgud(:)
call zfftifc(3,ngridg,-1,zfft1)
! (grad dcdgud).(grad rhodn) and (grad dcdgud).(grad rhoup)
do i=1,3
  zfft2(:)=0.d0
  do ig=1,ngvec
    ifg=igfft(ig)
    zfft2(ifg)=vgc(i,ig)*cmplx(-aimag(zfft1(ifg)),dble(zfft1(ifg)),8)
  end do
  call zfftifc(3,ngridg,1,zfft2)
  vcup(:)=vcup(:)-dble(zfft2(:))*gvdn(:,i)
  vcdn(:)=vcdn(:)-dble(zfft2(:))*gvup(:,i)
end do
deallocate(rfir,zfft1,zfft2)
return
end subroutine
!EOC