Merge lp:~nickpapior/siesta/4.1-phases into lp:siesta/4.1

Proposed by Nick Papior
Status: Merged
Merged at revision: 894
Proposed branch: lp:~nickpapior/siesta/4.1-phases
Merge into: lp:siesta/4.1
Diff against target: 548 lines (+97/-129) (has conflicts)
8 files modified
Src/pdos2g.F (+8/-15)
Src/pdos2k.F (+24/-30)
Src/pdos3g.F (+6/-14)
Src/pdos3k.F (+23/-29)
Src/pdosg.F (+2/-3)
Src/pdosk.F (+15/-17)
Src/pdoskp.F (+15/-21)
version.info (+4/-0)
Text conflict in version.info
To merge this branch: bzr merge lp:~nickpapior/siesta/4.1-phases
Reviewer Review Type Date Requested Status
Alberto Garcia Pending
Review via email: mp+343468@code.launchpad.net

Commit message

Clarified phase usages in pdos codes.

Now comments explain why the pdos*k routines used a different sign in the S(k) generation when performing the matrix multiplications.

I.e. this branch introduces nothing but clarified code and slight performance increase due to removal of a small array and lots of copying.

To post a comment you must log in.
lp:~nickpapior/siesta/4.1-phases updated
895. By Nick Papior

Amended pdoskp for the details outlined in the previous two commits

Now the phases are well explained in all pdos routines.

Revision history for this message
Alberto Garcia (albertog) wrote :

Good. I will merge.

Preview Diff

[H/L] Next/Prev Comment, [J/K] Next/Prev File, [N/P] Next/Prev Hunk
=== modified file 'Src/pdos2g.F'
--- Src/pdos2g.F 2017-09-27 13:56:46 +0000
+++ Src/pdos2g.F 2018-04-18 08:41:10 +0000
@@ -44,7 +44,7 @@
44C **** AUXILIARY *****************************************************44C **** AUXILIARY *****************************************************
45C real*8 haux(nuo,nuo) : Auxiliary space for the hamiltonian matrix45C real*8 haux(nuo,nuo) : Auxiliary space for the hamiltonian matrix
46C real*8 saux(nuo,nuo) : Auxiliary space for the overlap matrix46C real*8 saux(nuo,nuo) : Auxiliary space for the overlap matrix
47C real*8 psi(nuo,nuo) : Auxiliary space for the eigenvectors47C complex*16 psi(2,nuo,nuo) : Auxiliary space for the eigenvectors
48C **** OUTPUT ********************************************************48C **** OUTPUT ********************************************************
49C real*8 dtot(nhist,4) : Total density of states49C real*8 dtot(nhist,4) : Total density of states
50C real*8 dpr(nhist,nuo,4): Proyected density of states50C real*8 dpr(nhist,nuo,4): Proyected density of states
@@ -71,10 +71,11 @@
7171
72 real(dp)72 real(dp)
73 . H(maxnh,4), S(maxnh), E1, E2, sigma, eo(maxo*2),73 . H(maxnh,4), S(maxnh), E1, E2, sigma, eo(maxo*2),
74 . psi(2,2,nuotot,2*nuo), dtot(nhist,4), dpr(nhist,nuotot,4)74 . dtot(nhist,4), dpr(nhist,nuotot,4)
7575 complex(dp), target :: psi(2,nuotot,2*nuo)
76 complex(dp)76
77 . Haux(2,nuotot,2,nuo), Saux(2,nuotot,2,nuo), caux(2,nuotot)77 complex(dp) Haux(2,nuotot,2,nuo), Saux(2,nuotot,2,nuo)
78 complex(dp), pointer :: caux(:,:)
78 external cdiag79 external cdiag
7980
80C Internal variables ---------------------------------------------------81C Internal variables ---------------------------------------------------
@@ -181,11 +182,7 @@
181 diff = (ener - eo(ibandg))**2 / (sigma ** 2)182 diff = (ener - eo(ibandg))**2 / (sigma ** 2)
182 if (diff .gt. 15.0d0) cycle183 if (diff .gt. 15.0d0) cycle
183 gauss = exp(-diff)184 gauss = exp(-diff)
184 caux(:,:) = dcmplx(0.0_dp,0.0_dp)185 caux => psi(:,:,iband) ! c_{up,j}, c_{down,j}
185 do j=1,nuotot
186 caux(1,j) = dcmplx(psi(1,1,j,iband),psi(2,1,j,iband)) ! c_{j,up)
187 caux(2,j) = dcmplx(psi(1,2,j,iband),psi(2,2,j,iband)) ! c_{j,down)
188 enddo
189 do jo = 1, Bnuo186 do jo = 1, Bnuo
190 call LocalToGlobalOrb(jo,BNode,Nodes,juo)187 call LocalToGlobalOrb(jo,BNode,Nodes,juo)
191 do io = 1, nuotot188 do io = 1, nuotot
@@ -218,11 +215,7 @@
218 diff = (ener - eo(iband))**2 / (sigma ** 2)215 diff = (ener - eo(iband))**2 / (sigma ** 2)
219 if (diff .gt. 15.0d0) cycle216 if (diff .gt. 15.0d0) cycle
220 gauss = exp(-diff)217 gauss = exp(-diff)
221 caux(:,:) = dcmplx(0.0_dp,0.0_dp)218 caux => psi(:,:,iband) ! c_{up,j}, c_{down,j}
222 do j=1,nuotot
223 caux(1,j) = dcmplx(psi(1,1,j,iband),psi(2,1,j,iband)) ! c_{j,up)
224 caux(2,j) = dcmplx(psi(1,2,j,iband),psi(2,2,j,iband)) ! c_{j,down)
225 enddo
226 do jo = 1, nuotot219 do jo = 1, nuotot
227 do io = 1, nuotot220 do io = 1, nuotot
228 rpipj = dreal(caux(1,io)*dconjg(caux(1,jo)))*gauss221 rpipj = dreal(caux(1,io)*dconjg(caux(1,jo)))*gauss
229222
=== modified file 'Src/pdos2k.F'
--- Src/pdos2k.F 2017-09-27 13:56:46 +0000
+++ Src/pdos2k.F 2018-04-18 08:41:10 +0000
@@ -48,7 +48,7 @@
48C **** AUXILIARY *****************************************************48C **** AUXILIARY *****************************************************
49C REAL*8 HAUX(2,NUO,NUO) : Auxiliary space for the hamiltonian matrix49C REAL*8 HAUX(2,NUO,NUO) : Auxiliary space for the hamiltonian matrix
50C REAL*8 SAUX(2,NUO,NUO) : Auxiliary space for the overlap matrix50C REAL*8 SAUX(2,NUO,NUO) : Auxiliary space for the overlap matrix
51C REAL*8 PSI(2,NUO,NUO) : Auxiliary space for the eigenvectors51C COMPLEX*16 PSI(2,NUO,NUO) : Auxiliary space for the eigenvectors
52C **** OUTPUT ********************************************************52C **** OUTPUT ********************************************************
53C REAL*8 DTOT(NHIST,4) : Total density of states53C REAL*8 DTOT(NHIST,4) : Total density of states
54C REAL*8 DPR(NHIST,NUO,4): Proyected density of states54C REAL*8 DPR(NHIST,NUO,4): Proyected density of states
@@ -75,11 +75,11 @@
75 real(dp)75 real(dp)
76 . H(MAXNH,4), S(MAXNH), E1, E2, SIGMA,76 . H(MAXNH,4), S(MAXNH), E1, E2, SIGMA,
77 . XIJ(3,MAXNH), KPOINT(3,NK), EO(MAXO*2,NK),77 . XIJ(3,MAXNH), KPOINT(3,NK), EO(MAXO*2,NK),
78 . DTOT(NHIST,4), DPR(NHIST,NUOTOT,4), WK(NK),78 . DTOT(NHIST,4), DPR(NHIST,NUOTOT,4), WK(NK)
79 . psi(2,2,nuotot,2*nuo)79 complex(dp), target :: psi(2,nuotot,2*nuo)
8080
81 complex(dp) Haux(2,nuotot,2,nuo), Saux(2,nuotot,2,nuo), 81 complex(dp) Haux(2,nuotot,2,nuo), Saux(2,nuotot,2,nuo)
82 . caux(2,nuotot)82 complex(dp), pointer :: caux(:,:)
83 83
84 complex(dp) kphs84 complex(dp) kphs
8585
@@ -94,8 +94,7 @@
94 real(dp)94 real(dp)
95 . KXIJ, CKXIJ, SKXIJ, DELTA, ENER, DIFF, GAUSS, NORM, WKSUM95 . KXIJ, CKXIJ, SKXIJ, DELTA, ENER, DIFF, GAUSS, NORM, WKSUM
9696
97 complex(dp)97 complex(dp) :: D11, D12, D22
98 . D11, D12, D21, D22
9998
100 complex(dp), pointer :: Spr(:,:) => null()99 complex(dp), pointer :: Spr(:,:) => null()
101100
@@ -127,7 +126,7 @@
127 kxij = kpoint(1,ik) * xij(1,ind) +126 kxij = kpoint(1,ik) * xij(1,ind) +
128 . kpoint(2,ik) * xij(2,ind) +127 . kpoint(2,ik) * xij(2,ind) +
129 . kpoint(3,ik) * xij(3,ind)128 . kpoint(3,ik) * xij(3,ind)
130 kphs = cdexp(dcmplx(0.0_dp, -1.0_dp)*kxij)129 kphs = cdexp(dcmplx(0.0_dp, -kxij))
131130
132 Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind) * kphs131 Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind) * kphs
133 Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind) * kphs132 Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind) * kphs
@@ -159,7 +158,7 @@
159 kxij = kpoint(1,ik) * xij(1,ind) +158 kxij = kpoint(1,ik) * xij(1,ind) +
160 . kpoint(2,ik) * xij(2,ind) +159 . kpoint(2,ik) * xij(2,ind) +
161 . kpoint(3,ik) * xij(3,ind)160 . kpoint(3,ik) * xij(3,ind)
162 kphs = cdexp(dcmplx(0.0_dp, -1.0_dp)*kxij)161 kphs = cdexp(dcmplx(0.0_dp, -kxij))
163162
164 Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind) * kphs163 Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind) * kphs
165 Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind) * kphs164 Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind) * kphs
@@ -188,7 +187,10 @@
188 kxij = kpoint(1,ik) * xij(1,ind) +187 kxij = kpoint(1,ik) * xij(1,ind) +
189 . kpoint(2,ik) * xij(2,ind) +188 . kpoint(2,ik) * xij(2,ind) +
190 . kpoint(3,ik) * xij(3,ind)189 . kpoint(3,ik) * xij(3,ind)
191 kphs = cdexp(dcmplx(0.0_dp, 1.0_dp)*kxij)190! Since we are doing element wise multiplications (and not dot-products)
191! we might as well setup the transpose S(k)^T == S(-k) because this will
192! mean that we can do a simpler multiplication further down
193 kphs = cdexp(dcmplx(0.0_dp, kxij))
192 Spr(juo,iuo) = Spr(juo,iuo) + S(ind) * kphs194 Spr(juo,iuo) = Spr(juo,iuo) + S(ind) * kphs
193 enddo195 enddo
194 enddo196 enddo
@@ -225,23 +227,19 @@
225 ibandg = ibandg * 2 - mod(iband, 2) 227 ibandg = ibandg * 2 - mod(iband, 2)
226 diff = (ener - eo(ibandg,ik))**2 / (sigma ** 2)228 diff = (ener - eo(ibandg,ik))**2 / (sigma ** 2)
227 if (diff .gt. 15.0d0) cycle229 if (diff .gt. 15.0d0) cycle
228 gauss = exp(-diff)230 gauss = exp(-diff) * wk(ik)
229 caux(:,:) = dcmplx(0.0_dp,0.0_dp)231 caux => psi(:,:,iband) ! c_{j,up), c_{j,down}
230 do j=1,nuotot
231 caux(1,j) = dcmplx(psi(1,1,j,iband),psi(2,1,j,iband)) ! c_{j,up)
232 caux(2,j) = dcmplx(psi(1,2,j,iband),psi(2,2,j,iband)) ! c_{j,down)
233 enddo
234 do jo = 1, Bnuo232 do jo = 1, Bnuo
235 call LocalToGlobalOrb(jo,BNode,Nodes,juo)233 call LocalToGlobalOrb(jo,BNode,Nodes,juo)
236 do io = 1, nuotot234 do io = 1, nuotot
237 D11 = caux(1,io) * dconjg(caux(1,juo)) * Sloc(io,jo)235 D11 = caux(1,io) * dconjg(caux(1,juo)) * Sloc(io,jo)
238 D22 = caux(2,io) * dconjg(caux(2,juo)) * Sloc(io,jo)236 D22 = caux(2,io) * dconjg(caux(2,juo)) * Sloc(io,jo)
239 D12 = caux(1,io) * dconjg(caux(2,juo)) * Sloc(io,jo)237 D12 = caux(1,io) * dconjg(caux(2,juo)) * Sloc(io,jo)
240 D21 = caux(2,io) * dconjg(caux(1,juo)) * Sloc(io,jo)238 !D21 = caux(2,io) * dconjg(caux(1,juo)) * Sloc(io,jo)
241239
242 D11 = gauss*wk(ik)*D11240 D11 = gauss*D11
243 D22 = gauss*wk(ik)*D22241 D22 = gauss*D22
244 D12 = gauss*wk(ik)*D12242 D12 = gauss*D12
245243
246 dpr(ihist,juo,1) = dpr(ihist,juo,1) + dreal(D11)244 dpr(ihist,juo,1) = dpr(ihist,juo,1) + dreal(D11)
247 dpr(ihist,juo,2) = dpr(ihist,juo,2) + dreal(D22)245 dpr(ihist,juo,2) = dpr(ihist,juo,2) + dreal(D22)
@@ -265,22 +263,18 @@
265 do iband = 1, nuo*2263 do iband = 1, nuo*2
266 diff = (ener - eo(iband,ik))**2 / (sigma ** 2)264 diff = (ener - eo(iband,ik))**2 / (sigma ** 2)
267 if (diff .gt. 15.0d0) cycle265 if (diff .gt. 15.0d0) cycle
268 gauss = exp(-diff)266 gauss = exp(-diff) * wk(ik)
269 caux(:,:) = dcmplx(0.0_dp,0.0_dp)267 caux => psi(:,:,iband) ! c_{up,j), c_{down,j}
270 do j=1,nuotot
271 caux(1,j) = dcmplx(psi(1,1,j,iband),psi(2,1,j,iband)) ! c_{j,up)
272 caux(2,j) = dcmplx(psi(1,2,j,iband),psi(2,2,j,iband)) ! c_{j,down)
273 enddo
274 do io = 1, nuotot268 do io = 1, nuotot
275 do jo = 1, nuotot269 do jo = 1, nuotot
276 D11 = caux(1,io) * dconjg(caux(1,jo)) * Spr(io,jo)270 D11 = caux(1,io) * dconjg(caux(1,jo)) * Spr(io,jo)
277 D22 = caux(2,io) * dconjg(caux(2,jo)) * Spr(io,jo)271 D22 = caux(2,io) * dconjg(caux(2,jo)) * Spr(io,jo)
278 D12 = caux(1,io) * dconjg(caux(2,jo)) * Spr(io,jo)272 D12 = caux(1,io) * dconjg(caux(2,jo)) * Spr(io,jo)
279 D21 = caux(2,io) * dconjg(caux(1,jo)) * Spr(io,jo)273 !D21 = caux(2,io) * dconjg(caux(1,jo)) * Spr(io,jo)
280274
281 D11 = gauss*wk(ik)*D11275 D11 = gauss*D11
282 D22 = gauss*wk(ik)*D22276 D22 = gauss*D22
283 D12 = gauss*wk(ik)*D12277 D12 = gauss*D12
284278
285 dpr(ihist,jo,1) = dpr(ihist,jo,1) + dreal(D11)279 dpr(ihist,jo,1) = dpr(ihist,jo,1) + dreal(D11)
286 dpr(ihist,jo,2) = dpr(ihist,jo,2) + dreal(D22)280 dpr(ihist,jo,2) = dpr(ihist,jo,2) + dreal(D22)
287281
=== modified file 'Src/pdos3g.F'
--- Src/pdos3g.F 2017-09-27 13:56:46 +0000
+++ Src/pdos3g.F 2018-04-18 08:41:10 +0000
@@ -71,10 +71,10 @@
7171
72 real(dp)72 real(dp)
73 . H(maxnh,8), S(maxnh), E1, E2, sigma, eo(maxo*2),73 . H(maxnh,8), S(maxnh), E1, E2, sigma, eo(maxo*2),
74 . psi(2,2,nuotot,2*nuo), dtot(nhist,4), dpr(nhist,nuotot,4)74 . dtot(nhist,4), dpr(nhist,nuotot,4)
7575 complex(dp), target :: psi(2,nuotot,2*nuo)
76 complex(dp)76 complex(dp) Haux(2,nuotot,2,nuo), Saux(2,nuotot,2,nuo)
77 . Haux(2,nuotot,2,nuo), Saux(2,nuotot,2,nuo), caux(2,nuotot)77 complex(dp), pointer :: caux(:,:)
78 external cdiag78 external cdiag
7979
80C Internal variables ---------------------------------------------------80C Internal variables ---------------------------------------------------
@@ -180,11 +180,7 @@
180 diff = (ener - eo(ibandg))**2 / (sigma ** 2)180 diff = (ener - eo(ibandg))**2 / (sigma ** 2)
181 if (diff .gt. 15.0d0) cycle181 if (diff .gt. 15.0d0) cycle
182 gauss = exp(-diff)182 gauss = exp(-diff)
183 caux(:,:) = dcmplx(0.0_dp,0.0_dp)183 caux => psi(:,:,iband) ! c_{up,j}, c_{down,j}
184 do j=1,nuotot
185 caux(1,j) = dcmplx(psi(1,1,j,iband),psi(2,1,j,iband)) ! c_{j,up)
186 caux(2,j) = dcmplx(psi(1,2,j,iband),psi(2,2,j,iband)) ! c_{j,down)
187 enddo
188 do jo = 1, Bnuo184 do jo = 1, Bnuo
189 call LocalToGlobalOrb(jo,BNode,Nodes,juo)185 call LocalToGlobalOrb(jo,BNode,Nodes,juo)
190 do io = 1, nuotot186 do io = 1, nuotot
@@ -218,11 +214,7 @@
218 diff = (ener - eo(iband))**2 / (sigma ** 2)214 diff = (ener - eo(iband))**2 / (sigma ** 2)
219 if (diff .gt. 15.0d0) cycle215 if (diff .gt. 15.0d0) cycle
220 gauss = exp(-diff)216 gauss = exp(-diff)
221 caux(:,:) = dcmplx(0.0_dp,0.0_dp)217 caux => psi(:,:,iband) ! c_{up,j}, c_{down,j}
222 do j=1,nuotot
223 caux(1,j) = dcmplx(psi(1,1,j,iband),psi(2,1,j,iband)) ! c_{j,up)
224 caux(2,j) = dcmplx(psi(1,2,j,iband),psi(2,2,j,iband)) ! c_{j,down)
225 enddo
226 do jo = 1, nuotot218 do jo = 1, nuotot
227 do io = 1, nuotot219 do io = 1, nuotot
228 rpipj = dreal(caux(1,io)*dconjg(caux(1,jo)))*gauss220 rpipj = dreal(caux(1,io)*dconjg(caux(1,jo)))*gauss
229221
=== modified file 'Src/pdos3k.F'
--- Src/pdos3k.F 2017-09-27 13:56:46 +0000
+++ Src/pdos3k.F 2018-04-18 08:41:10 +0000
@@ -75,11 +75,11 @@
75 real(dp)75 real(dp)
76 . H(MAXNH,8), S(MAXNH), E1, E2, SIGMA,76 . H(MAXNH,8), S(MAXNH), E1, E2, SIGMA,
77 . XIJ(3,MAXNH), KPOINT(3,NK), EO(MAXO*2,NK),77 . XIJ(3,MAXNH), KPOINT(3,NK), EO(MAXO*2,NK),
78 . DTOT(NHIST,4), DPR(NHIST,NUOTOT,4), WK(NK),78 . DTOT(NHIST,4), DPR(NHIST,NUOTOT,4), WK(NK)
79 . psi(2,2,nuotot,2*nuo)79 complex(dp), target :: psi(2,nuotot,2*nuo)
8080
81 complex(dp) Haux(2,nuotot,2,nuo), Saux(2,nuotot,2,nuo), 81 complex(dp) Haux(2,nuotot,2,nuo), Saux(2,nuotot,2,nuo)
82 . caux(2,nuotot)82 complex(dp), pointer :: caux(:,:)
83 83
84 complex(dp) kphs84 complex(dp) kphs
8585
@@ -94,8 +94,7 @@
94 real(dp)94 real(dp)
95 . KXIJ, CKXIJ, SKXIJ, DELTA, ENER, DIFF, GAUSS, NORM, WKSUM95 . KXIJ, CKXIJ, SKXIJ, DELTA, ENER, DIFF, GAUSS, NORM, WKSUM
9696
97 complex(dp)97 complex(dp) :: D11, D12, D22
98 . D11, D12, D21, D22
9998
100 complex(dp), pointer :: Spr(:,:) => null()99 complex(dp), pointer :: Spr(:,:) => null()
101100
@@ -127,7 +126,7 @@
127 kxij = kpoint(1,ik) * xij(1,ind) +126 kxij = kpoint(1,ik) * xij(1,ind) +
128 . kpoint(2,ik) * xij(2,ind) +127 . kpoint(2,ik) * xij(2,ind) +
129 . kpoint(3,ik) * xij(3,ind)128 . kpoint(3,ik) * xij(3,ind)
130 kphs = cdexp(dcmplx(0.0_dp, -1.0_dp)*kxij)129 kphs = cdexp(dcmplx(0.0_dp, -kxij))
131130
132 Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind) * kphs131 Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind) * kphs
133 Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind) * kphs132 Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind) * kphs
@@ -159,7 +158,7 @@
159 kxij = kpoint(1,ik) * xij(1,ind) +158 kxij = kpoint(1,ik) * xij(1,ind) +
160 . kpoint(2,ik) * xij(2,ind) +159 . kpoint(2,ik) * xij(2,ind) +
161 . kpoint(3,ik) * xij(3,ind)160 . kpoint(3,ik) * xij(3,ind)
162 kphs = cdexp(dcmplx(0.0_dp, -1.0_dp)*kxij)161 kphs = cdexp(dcmplx(0.0_dp, -kxij))
163162
164 Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind) * kphs163 Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind) * kphs
165 Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind) * kphs164 Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind) * kphs
@@ -188,7 +187,10 @@
188 kxij = kpoint(1,ik) * xij(1,ind) +187 kxij = kpoint(1,ik) * xij(1,ind) +
189 . kpoint(2,ik) * xij(2,ind) +188 . kpoint(2,ik) * xij(2,ind) +
190 . kpoint(3,ik) * xij(3,ind)189 . kpoint(3,ik) * xij(3,ind)
191 kphs = cdexp(dcmplx(0.0_dp, 1.0_dp)*kxij)190! Since we are doing element wise multiplications (and not dot-products)
191! we might as well setup the transpose S(k)^T == S(-k) because this will
192! mean that we can do a simpler multiplication further down
193 kphs = cdexp(dcmplx(0.0_dp, kxij))
192 Spr(juo,iuo) = Spr(juo,iuo) + S(ind) * kphs194 Spr(juo,iuo) = Spr(juo,iuo) + S(ind) * kphs
193 enddo195 enddo
194 enddo196 enddo
@@ -225,23 +227,19 @@
225 ibandg = ibandg * 2 - mod(iband, 2) 227 ibandg = ibandg * 2 - mod(iband, 2)
226 diff = (ener - eo(ibandg,ik))**2 / (sigma ** 2)228 diff = (ener - eo(ibandg,ik))**2 / (sigma ** 2)
227 if (diff .gt. 15.0d0) cycle229 if (diff .gt. 15.0d0) cycle
228 gauss = exp(-diff)230 gauss = exp(-diff) * wk(ik)
229 caux(:,:) = dcmplx(0.0_dp,0.0_dp)231 caux => psi(:,:,iband) ! c_{up,j}, c_{down,j}
230 do j=1,nuotot
231 caux(1,j) = dcmplx(psi(1,1,j,iband),psi(2,1,j,iband)) ! c_{j,up)
232 caux(2,j) = dcmplx(psi(1,2,j,iband),psi(2,2,j,iband)) ! c_{j,down)
233 enddo
234 do jo = 1, Bnuo232 do jo = 1, Bnuo
235 call LocalToGlobalOrb(jo,BNode,Nodes,juo)233 call LocalToGlobalOrb(jo,BNode,Nodes,juo)
236 do io = 1, nuotot234 do io = 1, nuotot
237 D11 = caux(1,io) * dconjg(caux(1,juo)) * Sloc(io,jo)235 D11 = caux(1,io) * dconjg(caux(1,juo)) * Sloc(io,jo)
238 D22 = caux(2,io) * dconjg(caux(2,juo)) * Sloc(io,jo)236 D22 = caux(2,io) * dconjg(caux(2,juo)) * Sloc(io,jo)
239 D12 = caux(1,io) * dconjg(caux(2,juo)) * Sloc(io,jo)237 D12 = caux(1,io) * dconjg(caux(2,juo)) * Sloc(io,jo)
240 D21 = caux(2,io) * dconjg(caux(1,juo)) * Sloc(io,jo)238 !D21 = caux(2,io) * dconjg(caux(1,juo)) * Sloc(io,jo)
241239
242 D11 = gauss*wk(ik)*D11240 D11 = gauss*D11
243 D22 = gauss*wk(ik)*D22241 D22 = gauss*D22
244 D12 = gauss*wk(ik)*D12242 D12 = gauss*D12
245243
246 dpr(ihist,juo,1) = dpr(ihist,juo,1) + dreal(D11)244 dpr(ihist,juo,1) = dpr(ihist,juo,1) + dreal(D11)
247 dpr(ihist,juo,2) = dpr(ihist,juo,2) + dreal(D22)245 dpr(ihist,juo,2) = dpr(ihist,juo,2) + dreal(D22)
@@ -265,22 +263,18 @@
265 do iband = 1, nuo*2263 do iband = 1, nuo*2
266 diff = (ener - eo(iband,ik))**2 / (sigma ** 2)264 diff = (ener - eo(iband,ik))**2 / (sigma ** 2)
267 if (diff .gt. 15.0d0) cycle265 if (diff .gt. 15.0d0) cycle
268 gauss = exp(-diff)266 gauss = exp(-diff) * wk(ik)
269 caux(:,:) = dcmplx(0.0_dp,0.0_dp)267 caux => psi(:,:,iband) ! c_{up,j}, c_{down,j}
270 do j=1,nuotot
271 caux(1,j) = dcmplx(psi(1,1,j,iband),psi(2,1,j,iband)) ! c_{j,up)
272 caux(2,j) = dcmplx(psi(1,2,j,iband),psi(2,2,j,iband)) ! c_{j,down)
273 enddo
274 do io = 1, nuotot268 do io = 1, nuotot
275 do jo = 1, nuotot269 do jo = 1, nuotot
276 D11 = caux(1,io) * dconjg(caux(1,jo)) * Spr(io,jo)270 D11 = caux(1,io) * dconjg(caux(1,jo)) * Spr(io,jo)
277 D22 = caux(2,io) * dconjg(caux(2,jo)) * Spr(io,jo)271 D22 = caux(2,io) * dconjg(caux(2,jo)) * Spr(io,jo)
278 D12 = caux(1,io) * dconjg(caux(2,jo)) * Spr(io,jo)272 D12 = caux(1,io) * dconjg(caux(2,jo)) * Spr(io,jo)
279 D21 = caux(2,io) * dconjg(caux(1,jo)) * Spr(io,jo)273 !D21 = caux(2,io) * dconjg(caux(1,jo)) * Spr(io,jo)
280274
281 D11 = gauss*wk(ik)*D11275 D11 = gauss*D11
282 D22 = gauss*wk(ik)*D22276 D22 = gauss*D22
283 D12 = gauss*wk(ik)*D12277 D12 = gauss*D12
284278
285 dpr(ihist,jo,1) = dpr(ihist,jo,1) + dreal(D11)279 dpr(ihist,jo,1) = dpr(ihist,jo,1) + dreal(D11)
286 dpr(ihist,jo,2) = dpr(ihist,jo,2) + dreal(D22)280 dpr(ihist,jo,2) = dpr(ihist,jo,2) + dreal(D22)
287281
=== modified file 'Src/pdosg.F'
--- Src/pdosg.F 2017-07-14 12:17:28 +0000
+++ Src/pdosg.F 2018-04-18 08:41:10 +0000
@@ -190,7 +190,7 @@
190 if (diff .gt. 15.0d0) then190 if (diff .gt. 15.0d0) then
191 cycle191 cycle
192 else192 else
193 gauss = ( exp(-diff) )193 gauss = exp(-diff)
194 if (Node.eq.BNode) then194 if (Node.eq.BNode) then
195C Only add once to dtot - not everytime loop over processors is executed195C Only add once to dtot - not everytime loop over processors is executed
196 dtot(ihist,ispin) = dtot(ihist,ispin) + gauss196 dtot(ihist,ispin) = dtot(ihist,ispin) + gauss
@@ -224,10 +224,9 @@
224 if (diff .gt. 15.0d0) then224 if (diff .gt. 15.0d0) then
225 cycle225 cycle
226 else226 else
227 gauss = ( exp(-diff) )227 gauss = exp(-diff)
228 dtot(ihist,ispin) = dtot(ihist,ispin) + gauss228 dtot(ihist,ispin) = dtot(ihist,ispin) + gauss
229 do iuo = 1, nuotot229 do iuo = 1, nuotot
230C Solo para los Juo que satisfagan el criterio del record...
231 do juo = 1, nuotot230 do juo = 1, nuotot
232 pipj1 = psi(iuo,iband) * psi(juo,iband)231 pipj1 = psi(iuo,iband) * psi(juo,iband)
233 dpr(ihist,juo,ispin) = dpr(ihist,juo,ispin) + 232 dpr(ihist,juo,ispin) = dpr(ihist,juo,ispin) +
234233
=== modified file 'Src/pdosk.F'
--- Src/pdosk.F 2016-06-04 20:06:11 +0000
+++ Src/pdosk.F 2018-04-18 08:41:10 +0000
@@ -87,7 +87,7 @@
8787
88 real(dp)88 real(dp)
89 . kxij, Ckxij, Skxij, delta, ener, diff, pipj1, pipj2, 89 . kxij, Ckxij, Skxij, delta, ener, diff, pipj1, pipj2,
90 . pipjS1, pipjS2, gauss, norm, wksum90 . pipjS1, gauss, norm, wksum
9191
92#ifdef MPI92#ifdef MPI
93 integer ::93 integer ::
@@ -197,10 +197,11 @@
197 . kpoint(3,IK) * xij(3,ind) 197 . kpoint(3,IK) * xij(3,ind)
198 ckxij = cos(kxij)198 ckxij = cos(kxij)
199 skxij = sin(kxij)199 skxij = sin(kxij)
200C Calculates the hamiltonian and the overlap in k space200! Since we are doing element wise multiplications (and not dot-products)
201C H(k) = Sum(R) exp(i*k*R) * H(R)201! we might as well setup the transpose S(k)^T == S(-k) because this will
202! mean that we can do a simpler multiplication further down
202 Saux(1,juo,iuo) = Saux(1,juo,iuo) + S(ind) * ckxij203 Saux(1,juo,iuo) = Saux(1,juo,iuo) + S(ind) * ckxij
203 Saux(2,juo,iuo) = Saux(2,juo,iuo) - S(ind) * skxij204 Saux(2,juo,iuo) = Saux(2,juo,iuo) + S(ind) * skxij
204 enddo205 enddo
205 enddo206 enddo
206207
@@ -236,23 +237,22 @@
236 if (diff .gt. 15.0D0) then237 if (diff .gt. 15.0D0) then
237 cycle238 cycle
238 else239 else
239 gauss = ( EXP(-diff) )240 gauss = exp(-diff) * wk(ik)
240 if (Node.eq.BNode) then241 if (Node.eq.BNode) then
241C Only add once to dtot - not everytime loop over processors is executed242C Only add once to dtot - not everytime loop over processors is executed
242 dtot(ihist,ispin) = dtot(ihist,ispin) + gauss*WK(IK)243 dtot(ihist,ispin) = dtot(ihist,ispin) + gauss
243 endif244 endif
244 do jo = 1, Bnuo245 do jo = 1, Bnuo
245 call LocalToGlobalOrb(jo,BNode,Nodes,juo)246 call LocalToGlobalOrb(jo,BNode,Nodes,juo)
246 do iuo = 1, nuotot247 do iuo = 1, nuotot
247C Solo para los Juo que satisfagan el criterio del record...248! This is: psi(iuo) * psi(juo)^*
248 pipj1 = psi(1,iuo,iband) * psi(1,juo,iband) +249 pipj1 = psi(1,iuo,iband) * psi(1,juo,iband) +
249 . psi(2,iuo,iband) * psi(2,juo,iband)250 . psi(2,iuo,iband) * psi(2,juo,iband)
250 pipj2 = psi(1,iuo,iband) * psi(2,juo,iband) -251 pipj2 = - psi(1,iuo,iband) * psi(2,juo,iband) +
251 . psi(2,iuo,iband) * psi(1,juo,iband)252 . psi(2,iuo,iband) * psi(1,juo,iband)
252 pipjS1= pipj1*Sloc(1,iuo,JO)-pipj2*Sloc(2,iuo,JO)253 pipjS1= pipj1*Sloc(1,iuo,JO)-pipj2*Sloc(2,iuo,JO)
253 pipjS2= pipj1*Sloc(2,iuo,JO)+pipj2*Sloc(1,iuo,JO)254 dpr(ihist,juo,ispin)= dpr(ihist,juo,ispin) +
254 dpr(ihist,juo,ispin)= dpr(ihist,juo,ispin) + 255 . pipjS1*gauss
255 . pipjS1*gauss*WK(IK)
256 enddo256 enddo
257 enddo257 enddo
258 endif258 endif
@@ -275,19 +275,17 @@
275 if (diff .gt. 15.0d0) then275 if (diff .gt. 15.0d0) then
276 cycle276 cycle
277 else277 else
278 gauss = ( EXP(-diff) )278 gauss = exp(-diff) * wk(ik)
279 dtot(ihist,ispin) = dtot(ihist,ispin) + gauss*WK(IK)279 dtot(ihist,ispin) = dtot(ihist,ispin) + gauss
280 do iuo = 1, nuotot280 do iuo = 1, nuotot
281C Solo para los Juo que satisfagan el criterio del record...
282 do juo = 1, nuotot281 do juo = 1, nuotot
283 pipj1 = psi(1,iuo,iband) * psi(1,juo,iband) +282 pipj1 = psi(1,iuo,iband) * psi(1,juo,iband) +
284 . psi(2,iuo,iband) * psi(2,juo,iband)283 . psi(2,iuo,iband) * psi(2,juo,iband)
285 pipj2 = psi(1,iuo,iband) * psi(2,juo,iband) -284 pipj2 = - psi(1,iuo,iband) * psi(2,juo,iband) +
286 . psi(2,iuo,iband) * psi(1,juo,iband)285 . psi(2,iuo,iband) * psi(1,juo,iband)
287 pipjS1= pipj1*Saux(1,iuo,juo)-pipj2*Saux(2,iuo,juo)286 pipjS1= pipj1*Saux(1,iuo,juo)-pipj2*Saux(2,iuo,juo)
288 pipjS2= pipj1*Saux(2,iuo,juo)+pipj2*Saux(1,iuo,juo)
289 dpr(ihist,juo,ispin)= dpr(ihist,juo,ispin) + 287 dpr(ihist,juo,ispin)= dpr(ihist,juo,ispin) +
290 . pipjS1*gauss*WK(IK)288 . pipjS1*gauss
291 enddo289 enddo
292 enddo290 enddo
293 endif291 endif
294292
=== modified file 'Src/pdoskp.F'
--- Src/pdoskp.F 2016-04-05 09:54:38 +0000
+++ Src/pdoskp.F 2018-04-18 08:41:10 +0000
@@ -94,9 +94,9 @@
9494
95 real(dp)95 real(dp)
96 . kxij, Ckxij, Skxij, delta, ener, diff, pipj1, pipj2, 96 . kxij, Ckxij, Skxij, delta, ener, diff, pipj1, pipj2,
97 . pipjS1, pipjS2, gauss, norm, wksum97 . pipjS1, gauss, norm, wksum
9898
99 real(dp), dimension(:), pointer :: Snew, Dloc99 real(dp), dimension(:), pointer :: Snew
100 real(dp), dimension(:,:), pointer :: Hnew, xijloc100 real(dp), dimension(:,:), pointer :: Hnew, xijloc
101101
102#ifdef MPI102#ifdef MPI
@@ -118,10 +118,6 @@
118 call re_alloc( listhptrg, 1, nuotot, name='listhptrg',118 call re_alloc( listhptrg, 1, nuotot, name='listhptrg',
119 & routine='pdoskp' )119 & routine='pdoskp' )
120120
121C Find maximum value in numh and create local storage
122 nullify( Dloc )
123 call re_alloc( Dloc, 1, nuotot, name='Dloc', routine='pdoskp' )
124
125C Globalise numh121C Globalise numh
126 do io = 1,nuotot122 do io = 1,nuotot
127 call WhichNodeOrb(io,Nodes,BNode)123 call WhichNodeOrb(io,Nodes,BNode)
@@ -290,16 +286,17 @@
290 jo = listhg(ind)286 jo = listhg(ind)
291 iuo = indxuo(io)287 iuo = indxuo(io)
292 juo = indxuo(jo)288 juo = indxuo(jo)
293C Calculates the phases k*r_ij289C Calculate the phases k*r_ij
294 kxij = kpoint(1,ik) * xijloc(1,ind) +290 kxij = kpoint(1,ik) * xijloc(1,ind) +
295 . kpoint(2,ik) * xijloc(2,ind) +291 . kpoint(2,ik) * xijloc(2,ind) +
296 . kpoint(3,ik) * xijloc(3,ind) 292 . kpoint(3,ik) * xijloc(3,ind)
297 ckxij = cos(kxij)293 ckxij = cos(kxij)
298 skxij = sin(kxij)294 skxij = sin(kxij)
299C Calculates the hamiltonian and the overlap in k space295! Since we are doing element wise multiplications (and not dot-products)
300C H(k) = Sum(R) exp(i*k*R) * H(R)296! we might as well setup the transpose S(k)^T == S(-k) because this will
297! mean that we can do a simpler multiplication further down
301 Saux(1,juo,iuo) = Saux(1,juo,iuo) + Snew(ind) * ckxij298 Saux(1,juo,iuo) = Saux(1,juo,iuo) + Snew(ind) * ckxij
302 Saux(2,juo,iuo) = Saux(2,juo,iuo) - Snew(ind) * skxij299 Saux(2,juo,iuo) = Saux(2,juo,iuo) + Snew(ind) * skxij
303 enddo300 enddo
304 enddo301 enddo
305302
@@ -311,19 +308,17 @@
311 if (diff .gt. 15.0d0) then308 if (diff .gt. 15.0d0) then
312 cycle309 cycle
313 else310 else
314 gauss = ( exp(-diff) )311 gauss = exp(-diff) * wk(ik)
315 dtot(ihist,is) = dtot(ihist,is) + gauss*wk(ik)312 dtot(ihist,is) = dtot(ihist,is) + gauss
316 do iuo = 1, nuotot313 do juo = 1, nuotot
317C Solo para los Juo que satisfagan el criterio del record...314 do iuo = 1, nuotot
318 do juo = 1, nuotot315! This is: psi(iuo) * psi(juo)^*
319 pipj1 = psi(1,iuo,iband) * psi(1,juo,iband) +316 pipj1 = psi(1,iuo,iband) * psi(1,juo,iband) +
320 . psi(2,iuo,iband) * psi(2,juo,iband)317 . psi(2,iuo,iband) * psi(2,juo,iband)
321 pipj2 = psi(1,iuo,iband) * psi(2,juo,iband) -318 pipj2 = - psi(1,iuo,iband) * psi(2,juo,iband) +
322 . psi(2,iuo,iband) * psi(1,juo,iband)319 . psi(2,iuo,iband) * psi(1,juo,iband)
323 pipjS1= pipj1*Saux(1,iuo,juo)-pipj2*Saux(2,iuo,juo)320 pipjS1 = pipj1*Saux(1,iuo,juo)-pipj2*Saux(2,iuo,juo)
324 pipjS2= pipj1*Saux(2,iuo,juo)+pipj2*Saux(1,iuo,juo)321 dpr(ihist,juo,is) = dpr(ihist,juo,is) + pipjS1*gauss
325 dpr(ihist,juo,is) = dpr(ihist,juo,is) +
326 . pipjS1*gauss*wk(ik)
327 enddo322 enddo
328 enddo323 enddo
329 endif324 endif
@@ -339,7 +334,6 @@
339 call de_alloc( xijloc, name='xijloc' )334 call de_alloc( xijloc, name='xijloc' )
340 call de_alloc( Hnew, name='Hnew' )335 call de_alloc( Hnew, name='Hnew' )
341 call de_alloc( Snew, name='Snew' )336 call de_alloc( Snew, name='Snew' )
342 call de_alloc( Dloc, name='Dloc' )
343 call de_alloc( listhg, name='listhg' )337 call de_alloc( listhg, name='listhg' )
344 call de_alloc( listhptrg, name='listhptrg' )338 call de_alloc( listhptrg, name='listhptrg' )
345 call de_alloc( numhg, name='numhg' )339 call de_alloc( numhg, name='numhg' )
346340
=== modified file 'version.info'
--- version.info 2018-04-17 13:06:00 +0000
+++ version.info 2018-04-18 08:41:10 +0000
@@ -1,1 +1,5 @@
1<<<<<<< TREE
1siesta-4.1--8932siesta-4.1--893
3=======
4siesta-4.1--892--phase-3
5>>>>>>> MERGE-SOURCE

Subscribers

People subscribed via source and target branches