Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 : MODULE m_vgen_xcpot
7 :
8 : USE m_juDFT
9 :
10 : CONTAINS
11 :
12 340 : SUBROUTINE vgen_xcpot(hybrid, input, xcpot, dimension, atoms, sphhar, stars, vacuum, sym, &
13 : obsolete, cell, oneD, sliceplot, mpi, noco, den, denRot, EnergyDen, vTot, vx, results)
14 :
15 : ! ***********************************************************
16 : ! FLAPW potential generator *
17 : ! ***********************************************************
18 : ! calculates the density-potential integrals needed for the
19 : ! total energy
20 : ! TE_VCOUL : charge density-coulomb potential integral
21 : ! TE_VEFF: charge density-effective potential integral
22 : ! TE_EXC : charge density-ex-corr.energy density integral
23 : ! ***********************************************************
24 :
25 : USE m_types
26 : USE m_constants
27 : USE m_intnv
28 : USE m_vmt_xc
29 : USE m_vvacxc
30 : USE m_vvacxcg
31 : USE m_vis_xc
32 : USE m_checkdopall
33 : USE m_cdn_io
34 : USE m_convol
35 : use m_cdntot
36 :
37 : IMPLICIT NONE
38 :
39 : CLASS(t_xcpot), INTENT(INOUT) :: xcpot
40 : TYPE(t_hybrid), INTENT(IN) :: hybrid
41 : TYPE(t_mpi), INTENT(IN) :: mpi
42 : TYPE(t_dimension), INTENT(IN) :: dimension
43 : TYPE(t_oneD), INTENT(IN) :: oneD
44 : TYPE(t_obsolete), INTENT(IN) :: obsolete
45 : TYPE(t_sliceplot), INTENT(IN) :: sliceplot
46 : TYPE(t_input), INTENT(IN) :: input
47 : TYPE(t_vacuum), INTENT(IN) :: vacuum
48 : TYPE(t_noco), INTENT(IN) :: noco
49 : TYPE(t_sym), INTENT(IN) :: sym
50 : TYPE(t_stars), INTENT(IN) :: stars
51 : TYPE(t_cell), INTENT(IN) :: cell
52 : TYPE(t_sphhar), INTENT(IN) :: sphhar
53 : TYPE(t_atoms), INTENT(IN) :: atoms
54 : TYPE(t_potden), INTENT(IN) :: den, denRot, EnergyDen
55 : TYPE(t_potden), INTENT(INOUT) :: vTot, vx
56 : TYPE(t_results), INTENT(INOUT), OPTIONAL :: results
57 :
58 : ! Local type instances
59 340 : TYPE(t_potden) :: workDen, exc, veff
60 : real, allocatable :: tmp_mt(:,:,:), tmp_is(:,:)
61 : ! Local Scalars
62 : INTEGER ifftd, ifftd2, ifftxc3d, ispin, i
63 : #ifdef CPP_MPI
64 : include 'mpif.h'
65 : integer:: ierr
66 : #endif
67 :
68 340 : CALL exc%init(stars, atoms, sphhar, vacuum, noco, 1, 1) !one spin only
69 340 : ALLOCATE (exc%pw_w(stars%ng3, 1)); exc%pw_w = 0.0
70 340 : IF (PRESENT(results)) THEN
71 340 : CALL veff%init(stars, atoms, sphhar, vacuum, noco, input%jspins, 1)
72 : #ifndef CPP_OLDINTEL
73 340 : ALLOCATE (veff%pw_w, mold=veff%pw)
74 : #else
75 : ALLOCATE (veff%pw_w(size(veff%pw, 1), size(veff%pw, 2)))
76 : #endif
77 : ENDIF
78 :
79 : ! exchange correlation potential
80 :
81 : ! vacuum region
82 340 : IF (mpi%irank == 0) THEN
83 170 : IF (input%film) THEN
84 5 : CALL timestart("Vxc in vacuum")
85 :
86 5 : ifftd2 = 9*stars%mx1*stars%mx2
87 5 : IF (oneD%odi%d1) ifftd2 = 9*stars%mx3*oneD%odi%M
88 :
89 5 : IF (.NOT. xcpot%needs_grad()) THEN ! LDA
90 :
91 0 : IF (.NOT. oneD%odi%d1) THEN
92 0 : CALL vvacxc(ifftd2, stars, vacuum, xcpot, input, noco, Den, vTot, exc)
93 : ELSE
94 0 : CALL judft_error("OneD broken")
95 : ! CALL vvacxc(stars,oneD%M,vacuum,odi%n2d,dimension,ifftd2,&
96 : ! xcpot,input,odi%nq2,odi%nst2,den,noco,odi%kimax2%igf,&
97 : ! odl%pgf,vTot%vacxy,vTot%vacz,excxy,excz)
98 : END IF
99 : ELSE ! GGA
100 5 : IF (oneD%odi%d1) THEN
101 0 : CALL judft_error("OneD broken")
102 : ! CALL vvacxcg(ifftd2,stars,vacuum,noco,oneD,&
103 : ! cell,xcpot,input,obsolete,workDen, ichsmrg,&
104 : ! vTot%vacxy,vTot%vacz,rhmn, exc%vacxy,exc%vacz)
105 :
106 : ELSE
107 5 : CALL vvacxcg(ifftd2, stars, vacuum, noco, oneD, cell, xcpot, input, obsolete, Den, vTot, exc)
108 : END IF
109 : END IF
110 5 : CALL timestop("Vxc in vacuum")
111 : END IF
112 :
113 : ! interstitial region
114 170 : CALL timestart("Vxc in interstitial")
115 :
116 170 : IF ((.NOT. obsolete%lwb) .OR. (.not. xcpot%needs_grad())) THEN
117 : ! no White-Bird-trick
118 170 : CALL vis_xc(stars, sym, cell, den, xcpot, input, noco, EnergyDen, vTot, vx, exc)
119 :
120 : ELSE
121 : ! White-Bird-trick
122 0 : WRITE (6, '(a)') "W+B trick cancelled out. visxcwb uses at present common block cpgft3.", &
123 0 : "visxcwb needs to be reprogrammed according to visxcg.f"
124 0 : CALL juDFT_error("visxcwb", calledby="vgen")
125 : END IF
126 :
127 170 : CALL timestop("Vxc in interstitial")
128 : END IF !irank==0
129 :
130 : !
131 : ! ------------------------------------------
132 : ! ----> muffin tin spheres region
133 :
134 340 : IF (mpi%irank == 0) THEN
135 170 : CALL timestart("Vxc in MT")
136 : END IF
137 :
138 : CALL vmt_xc(mpi, sphhar, atoms, den, xcpot, input, sym, &
139 340 : EnergyDen, vTot, vx, exc)
140 :
141 : ! add MT EXX potential to vr
142 340 : IF (mpi%irank == 0) THEN
143 170 : CALL timestop("Vxc in MT")
144 :
145 : ! check continuity of total potential
146 170 : IF (input%vchk) CALL checkDOPAll(input, dimension, sphhar, stars, atoms, sym, vacuum, oneD, cell, vTot, 1)
147 :
148 : ! TOTAL
149 170 : IF (PRESENT(results)) THEN
150 : ! CALCULATE THE INTEGRAL OF n1*Veff1 + n2*Veff2
151 : ! Veff = Vcoulomb + Vxc
152 170 : IF (noco%l_noco) THEN
153 124 : workDen = denRot
154 : ELSE
155 46 : workden = den
156 : END IF
157 :
158 170 : veff = vTot
159 170 : IF (xcpot%is_hybrid() .AND. hybrid%l_subvxc) THEN
160 0 : DO ispin = 1, input%jspins
161 0 : CALL convol(stars, vx%pw_w(:, ispin), vx%pw(:, ispin), stars%ufft)
162 : END DO
163 0 : veff%pw = vTot%pw - xcpot%get_exchange_weight()*vx%pw
164 0 : veff%pw_w = vTot%pw_w - xcpot%get_exchange_weight()*vx%pw_w
165 0 : veff%mt = vTot%mt - xcpot%get_exchange_weight()*vx%mt
166 0 : exc%pw = exc%pw - xcpot%get_exchange_weight()*exc%pw
167 0 : exc%pw_w = exc%pw_w - xcpot%get_exchange_weight()*exc%pw_w
168 0 : exc%mt = exc%mt - xcpot%get_exchange_weight()*exc%mt
169 : END IF
170 :
171 474 : DO ispin = 1, input%jspins
172 303114 : DO i = 1, stars%ng3
173 302640 : vx%pw(i, ispin) = vx%pw(i, ispin)/stars%nstr(i)
174 302944 : vx%pw_w(i, ispin) = vx%pw_w(i, ispin)/stars%nstr(i)
175 : END DO
176 : END DO
177 :
178 170 : results%te_veff = 0.0
179 474 : DO ispin = 1, input%jspins
180 304 : WRITE (6, FMT=8050) ispin
181 : 8050 FORMAT(/, 10x, 'density-effective potential integrals for spin ', i2,/)
182 474 : CALL int_nv(ispin, stars, vacuum, atoms, sphhar, cell, sym, input, oneD, veff, workden, results%te_veff)
183 : END DO
184 :
185 170 : WRITE (6, FMT=8060) results%te_veff
186 : 8060 FORMAT(/, 10x, 'total density-effective potential integral :', t40, ES20.10)
187 :
188 : ! CALCULATE THE INTEGRAL OF n*exc
189 :
190 : ! perform spin summation of charge densities for the calculation of Exc
191 170 : CALL workden%sum_both_spin()
192 :
193 170 : WRITE (6, FMT=8070)
194 : 8070 FORMAT(/, 10x, 'charge density-energy density integrals',/)
195 :
196 170 : results%te_exc = 0.0
197 170 : CALL int_nv(1, stars, vacuum, atoms, sphhar, cell, sym, input, oneD, exc, workDen, results%te_exc)
198 170 : WRITE (6, FMT=8080) results%te_exc
199 :
200 : 8080 FORMAT(/, 10x, 'total charge density-energy density integral :', t40, ES20.10)
201 : END IF
202 : END IF ! mpi%irank == 0
203 :
204 340 : END SUBROUTINE vgen_xcpot
205 :
206 : END MODULE m_vgen_xcpot
|
