Line data Source code
1 : MODULE m_greensfCalcImagPart
2 :
3 : USE m_types
4 : USE m_constants
5 : USE m_tetrahedronInit
6 : USE m_juDFT
7 :
8 : IMPLICIT NONE
9 :
10 : CONTAINS
11 :
12 456 : SUBROUTINE greensfCalcImagPart_single_kpt(ikpt,ikpt_i,ev_list,spin_ind,gfinp,atoms,input,kpts,noco,mpi,&
13 : results,greensfBZintCoeffs,greensfImagPart)
14 :
15 : INTEGER, INTENT(IN) :: ikpt, ikpt_i
16 : INTEGER, intent(IN) :: ev_list(:)
17 : INTEGER, INTENT(IN) :: spin_ind
18 : TYPE(t_gfinp), INTENT(IN) :: gfinp
19 : TYPE(t_atoms), INTENT(IN) :: atoms
20 : TYPE(t_input), INTENT(IN) :: input
21 : TYPE(t_kpts), INTENT(IN) :: kpts
22 : TYPE(t_noco), INTENT(IN) :: noco
23 : TYPE(t_mpi), INTENT(IN) :: mpi
24 : TYPE(t_results), INTENT(IN) :: results
25 : TYPE(t_greensfBZintCoeffs),INTENT(IN) :: greensfBZintCoeffs
26 : TYPE(t_greensfImagPart), INTENT(INOUT) :: greensfImagPart
27 :
28 : INTEGER :: nBands,jsp,i_gf,nLO,imatSize
29 : INTEGER :: l,lp,m,mp,iBand,ie,j,eGrid_start,eGrid_end
30 : INTEGER :: indUnique,i_elem,imat,iLO,iLOp,i_elemLO,i_elem_imag,i_elemLO_imag
31 : LOGICAL :: l_zero,l_sphavg,l_kresolved_int,l_kresolved
32 : REAL :: del,eb,wtkpt
33 : COMPLEX :: fac,weight
34 456 : REAL, ALLOCATABLE :: eig(:)
35 456 : REAL, ALLOCATABLE :: eMesh(:)
36 456 : COMPLEX, ALLOCATABLE :: imag(:,:)
37 456 : REAL, ALLOCATABLE :: weights(:,:)
38 456 : INTEGER, ALLOCATABLE :: indBound(:,:)
39 :
40 : !Get the information on the real axis energy mesh
41 456 : CALL gfinp%eMesh(results%ef,del,eb,eMesh=eMesh)
42 :
43 : !Spin degeneracy factors
44 456 : fac = 2.0/input%jspins
45 :
46 456 : nBands = SIZE(ev_list)
47 456 : jsp = MERGE(1,spin_ind,noco%l_noco)
48 24608 : eig = results%eig(ev_list,ikpt,jsp)
49 :
50 768 : SELECT CASE(input%bz_integration)
51 : CASE(BZINT_METHOD_HIST) !Histogram method
52 312 : wtkpt = kpts%wtkpt(ikpt)
53 : CASE(BZINT_METHOD_TETRA) !Tetrahedron method
54 144 : CALL timestart("Green's Function: TetrahedronWeights")
55 34524968 : ALLOCATE(weights(SIZE(eMesh),nBands),source=0.0)
56 15904 : ALLOCATE(indBound(nBands,2),source=0)
57 : CALL tetrahedronInit(kpts,input,ikpt,results%eig(ev_list,:,jsp),nBands,eMesh,&
58 67808 : weights,bounds=indBound,dos=.TRUE.)
59 34524536 : weights = - pi_const * fac * weights
60 144 : CALL timestop("Green's Function: TetrahedronWeights")
61 : CASE DEFAULT
62 : CALL juDFT_error("Invalid Brillouin-zone integration mode for Green's Functions",&
63 456 : hint="Choose either hist or tetra",calledby="greensfCalcImagPart")
64 : END SELECT
65 :
66 456 : CALL timestart("Green's Function: Imaginary Part")
67 : !Loop over Green's Function elements
68 10660 : DO i_gf = 1, gfinp%n
69 :
70 : !Get the information about the current element
71 10204 : l = gfinp%elem(i_gf)%l
72 10204 : lp = gfinp%elem(i_gf)%lp
73 10204 : l_sphavg = gfinp%elem(i_gf)%l_sphavg
74 10204 : l_kresolved_int = gfinp%elem(i_gf)%l_kresolved_int
75 10204 : l_kresolved = gfinp%elem(i_gf)%l_kresolved
76 :
77 :
78 10204 : IF(.NOT.gfinp%isUnique(i_gf, distinct_kresolved_int=.TRUE.)) CYCLE
79 :
80 5160 : i_elem = gfinp%uniqueElements(atoms,max_index=i_gf,l_sphavg=l_sphavg)
81 5160 : i_elemLO = gfinp%uniqueElements(atoms,max_index=i_gf,lo=.TRUE.,l_sphavg=l_sphavg)
82 :
83 5160 : i_elem_imag = gfinp%uniqueElements(atoms,max_index=i_gf,l_sphavg=l_sphavg,l_kresolved_int=l_kresolved_int)
84 5160 : i_elemLO_imag = gfinp%uniqueElements(atoms,max_index=i_gf,lo=.TRUE.,l_sphavg=l_sphavg,l_kresolved_int=l_kresolved_int)
85 :
86 5160 : nLO = 0
87 5160 : imatSize = 1
88 5160 : IF(.NOT.l_sphavg) THEN
89 120 : imatSize = 4
90 120 : nLO = gfinp%elem(i_gf)%countLOs(atoms)
91 120 : IF(nLO/=0) THEN
92 80 : imatSize = 4+4*nLO+nLO**2
93 : ENDIF
94 : ENDIF
95 :
96 : !$OMP parallel default(none) &
97 : !$OMP shared(input,gfinp,greensfBZintCoeffs,greensfImagPart) &
98 : !$OMP shared(i_elem,i_elemLO,i_elem_imag,i_elemLO_imag,nLO,l,lp,nBands,eMesh,l_sphavg,imatSize,l_kresolved_int,l_kresolved)&
99 : !$OMP shared(del,eb,eig,weights,indBound,fac,wtkpt,spin_ind,ikpt_i) &
100 5616 : !$OMP private(ie,iLO,iLOp,imat,m,mp,iBand,j,eGrid_start,eGrid_end,weight,imag,l_zero)
101 : ALLOCATE(imag(SIZE(eMesh),imatSize),source=cmplx_0)
102 : !$OMP do collapse(2)
103 : DO mp = -lp, lp
104 : DO m = -l, l
105 : imag = cmplx_0
106 : DO iBand = 1, nBands
107 :
108 : !Check for a non-zero weight in the energy interval
109 : l_zero = .TRUE.
110 : SELECT CASE(input%bz_integration)
111 : CASE(BZINT_METHOD_HIST) !Histogram Method
112 : j = FLOOR((eig(iBand)-eb)/del)+1
113 : IF(j.LE.SIZE(eMesh).AND.j.GE.1) l_zero = .FALSE.
114 : eGrid_start = j
115 : eGrid_end = j
116 : CASE(BZINT_METHOD_TETRA) !Tetrahedron method
117 : IF(ANY(ABS(weights(indBound(iBand,1):indBound(iBand,2),iBand)).GT.1e-14)) l_zero = .FALSE.
118 : eGrid_start = indBound(iBand,1)
119 : eGrid_end = indBound(iBand,2)
120 : CASE DEFAULT
121 : END SELECT
122 :
123 : IF(l_zero) CYCLE !No non-zero weight for this band
124 :
125 : IF(eGrid_start==eGrid_end) THEN
126 : ie = eGrid_start
127 :
128 : SELECT CASE(input%bz_integration)
129 : CASE(BZINT_METHOD_HIST) !Histogram Method
130 : IF(l_kresolved) THEN
131 : weight = -fac * pi_const/del
132 : ELSE
133 : weight = -fac * pi_const * wtkpt/del
134 : ENDIF
135 : CASE(BZINT_METHOD_TETRA) !Tetrahedron method
136 : weight = weights(ie,iBand)
137 : CASE DEFAULT
138 : END SELECT
139 :
140 : IF(l_sphavg) THEN
141 : imag(ie,1) = imag(ie,1) + weight * greensfBZintCoeffs%sphavg(iBand,m,mp,i_elem)
142 : ELSE
143 : imag(ie,1) = imag(ie,1) + weight * greensfBZintCoeffs%uu(iBand,m,mp,i_elem)
144 : imag(ie,2) = imag(ie,2) + weight * greensfBZintCoeffs%dd(iBand,m,mp,i_elem)
145 : imag(ie,3) = imag(ie,3) + weight * greensfBZintCoeffs%ud(iBand,m,mp,i_elem)
146 : imag(ie,4) = imag(ie,4) + weight * greensfBZintCoeffs%du(iBand,m,mp,i_elem)
147 : IF(nLO>0) THEN
148 : imat = 0
149 : DO iLO = 1, nLO
150 : imat = imat + 4
151 : imag(ie,imat+1) = imag(ie,imat+1) + weight * greensfBZintCoeffs%uulo(iBand,m,mp,iLO,i_elemLO)
152 : imag(ie,imat+2) = imag(ie,imat+2) + weight * greensfBZintCoeffs%ulou(iBand,m,mp,iLO,i_elemLO)
153 : imag(ie,imat+3) = imag(ie,imat+3) + weight * greensfBZintCoeffs%dulo(iBand,m,mp,iLO,i_elemLO)
154 : imag(ie,imat+4) = imag(ie,imat+4) + weight * greensfBZintCoeffs%ulod(iBand,m,mp,iLO,i_elemLO)
155 : ENDDO
156 : imat = 0
157 : DO iLO = 1, nLO
158 : DO iLOp = 1, nLO
159 : imat = imat + 1
160 : imag(ie,4 + 4*nLO+imat) = imag(ie,4 + 4*nLO+imat) + weight * greensfBZintCoeffs%uloulop(iBand,m,mp,imat,i_elemLO)
161 : ENDDO
162 : ENDDO
163 : ENDIF
164 : ENDIF
165 : ELSE
166 : IF(l_sphavg) THEN
167 : imag(eGrid_start:eGrid_end,1) = imag(eGrid_start:eGrid_end,1) + weights(eGrid_start:eGrid_end,iBand)&
168 : * greensfBZintCoeffs%sphavg(iBand,m,mp,i_elem)
169 : ELSE
170 : imag(eGrid_start:eGrid_end,1) = imag(eGrid_start:eGrid_end,1) + weights(eGrid_start:eGrid_end,iBand)&
171 : * greensfBZintCoeffs%uu(iBand,m,mp,i_elem)
172 : imag(eGrid_start:eGrid_end,2) = imag(eGrid_start:eGrid_end,2) + weights(eGrid_start:eGrid_end,iBand)&
173 : * greensfBZintCoeffs%dd(iBand,m,mp,i_elem)
174 : imag(eGrid_start:eGrid_end,3) = imag(eGrid_start:eGrid_end,3) + weights(eGrid_start:eGrid_end,iBand)&
175 : * greensfBZintCoeffs%ud(iBand,m,mp,i_elem)
176 : imag(eGrid_start:eGrid_end,4) = imag(eGrid_start:eGrid_end,4) + weights(eGrid_start:eGrid_end,iBand)&
177 : * greensfBZintCoeffs%du(iBand,m,mp,i_elem)
178 : IF(nLO>0) THEN
179 : imat = 0
180 : DO iLO = 1, nLO
181 : imat = imat + 4
182 : imag(eGrid_start:eGrid_end,imat+1) = imag(eGrid_start:eGrid_end,imat+1) + weights(eGrid_start:eGrid_end,iBand) &
183 : * greensfBZintCoeffs%uulo(iBand,m,mp,iLO,i_elemLO)
184 : imag(eGrid_start:eGrid_end,imat+2) = imag(eGrid_start:eGrid_end,imat+2) + weights(eGrid_start:eGrid_end,iBand) &
185 : * greensfBZintCoeffs%ulou(iBand,m,mp,iLO,i_elemLO)
186 : imag(eGrid_start:eGrid_end,imat+3) = imag(eGrid_start:eGrid_end,imat+3) + weights(eGrid_start:eGrid_end,iBand) &
187 : * greensfBZintCoeffs%dulo(iBand,m,mp,iLO,i_elemLO)
188 : imag(eGrid_start:eGrid_end,imat+4) = imag(eGrid_start:eGrid_end,imat+4) + weights(eGrid_start:eGrid_end,iBand) &
189 : * greensfBZintCoeffs%ulod(iBand,m,mp,iLO,i_elemLO)
190 : ENDDO
191 : imat = 0
192 : DO iLO = 1, nLO
193 : DO iLOp = 1, nLO
194 : imat = imat + 1
195 : imag(eGrid_start:eGrid_end,4 + 4*nLO+imat) = imag(eGrid_start:eGrid_end,4 + 4*nLO+imat) + weights(eGrid_start:eGrid_end,iBand) &
196 : * greensfBZintCoeffs%uloulop(iBand,m,mp,imat,i_elemLO)
197 : ENDDO
198 : ENDDO
199 : ENDIF
200 : ENDIF
201 : ENDIF
202 :
203 : ENDDO!ib
204 : IF(l_sphavg.AND..NOT.l_kresolved_int) THEN
205 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,1),1,greensfImagPart%sphavg(:,m,mp,i_elem_imag,spin_ind),1)
206 : ELSE IF(.NOT.l_kresolved_int) THEN
207 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,1),1,greensfImagPart%uu(:,m,mp,i_elem_imag,spin_ind),1)
208 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,2),1,greensfImagPart%dd(:,m,mp,i_elem_imag,spin_ind),1)
209 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,3),1,greensfImagPart%ud(:,m,mp,i_elem_imag,spin_ind),1)
210 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,4),1,greensfImagPart%du(:,m,mp,i_elem_imag,spin_ind),1)
211 :
212 : IF(nLO>0) THEN
213 : imat = 0
214 : DO iLO = 1, nLO
215 : imat = imat + 4
216 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,imat+1),1,greensfImagPart%uulo(:,m,mp,iLO,i_elemLO_imag,spin_ind),1)
217 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,imat+2),1,greensfImagPart%ulou(:,m,mp,iLO,i_elemLO_imag,spin_ind),1)
218 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,imat+3),1,greensfImagPart%dulo(:,m,mp,iLO,i_elemLO_imag,spin_ind),1)
219 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,imat+4),1,greensfImagPart%ulod(:,m,mp,iLO,i_elemLO_imag,spin_ind),1)
220 : ENDDO
221 : imat = 0
222 : DO iLO = 1, nLO
223 : DO iLOp = 1, nLO
224 : imat = imat + 1
225 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,4 + 4*nLO+imat),1,greensfImagPart%uloulop(:,m,mp,iLO,iLOp,i_elemLO_imag,spin_ind),1)
226 : ENDDO
227 : ENDDO
228 : ENDIF
229 : ELSE
230 : CALL zaxpy(SIZE(eMesh),cmplx_1,imag(:,1),1,greensfImagPart%sphavg_k(:,m,mp,i_elem_imag,spin_ind,ikpt_i),1)
231 : ENDIF
232 :
233 : ENDDO!m
234 : ENDDO!mp
235 : !$OMP end do
236 : DEALLOCATE(imag)
237 : !$OMP end parallel
238 : ENDDO!i_gf
239 456 : CALL timestop("Green's Function: Imaginary Part")
240 :
241 456 : IF(input%bz_integration==BZINT_METHOD_TETRA) DEALLOCATE(weights,indBound)
242 456 : end subroutine
243 :
244 : END MODULE m_greensfCalcImagPart
|
