• R/O
  • HTTP
  • SSH
  • HTTPS

libtetrabz: コミット

テトラへドロン法ライブラリ


コミットメタ情報

リビジョン4827df206f4f771e0aed915c5f87a244a7dde778 (tree)
日時2018-06-14 11:10:36
作者mitsuaki1987 <kawamitsuaki@gmai...>
コミッターmitsuaki1987

ログメッセージ

Use bi-linear interpolation; it may better for the extleamly coarse k-grid.

変更サマリ

差分

--- a/doc/index.html
+++ b/doc/index.html
@@ -8,7 +8,7 @@
88 <meta http-equiv="Pragma" content="no-cache">
99 <meta http-equiv="Cache-Control" content="no-cache">
1010
11- <title>Project Web of libtetrabz on OSDN Web space</title>
11+ <title>libTetraBZ Web page</title>
1212 </head>
1313
1414 <body bgcolor="CCFFCC">
--- a/doc/index.html.ja
+++ b/doc/index.html.ja
@@ -8,7 +8,7 @@
88 <meta http-equiv="Pragma" content="no-cache">
99 <meta http-equiv="Cache-Control" content="no-cache">
1010
11- <title>Project Web of libtetrabz on OSDN Web space</title>
11+ <title>libTetraBZ Web page</title>
1212 </head>
1313
1414 <body bgcolor="CCFFCC">
--- a/src/libtetrabz_common.F90
+++ b/src/libtetrabz_common.F90
@@ -336,103 +336,40 @@ END SUBROUTINE libtetrabz_sort
336336 !
337337 ! Linear interpolation
338338 !
339-SUBROUTINE libtetrabz_interpol_indx(nintp,ng,kvec,kintp,wintp)
339+SUBROUTINE libtetrabz_interpol_indx(ng,kvec,kintp,wintp)
340340 !
341341 IMPLICIT NONE
342342 !
343- INTEGER,INTENT(in) :: nintp, ng(3)
343+ INTEGER,INTENT(in) :: ng(3)
344344 REAL(8),INTENT(in) :: kvec(3)
345- INTEGER,INTENT(out) :: kintp(nintp)
346- REAL(8),INTENT(out) :: wintp(nintp)
345+ INTEGER,INTENT(out) :: kintp(8)
346+ REAL(8),INTENT(out) :: wintp(8)
347347 !
348- INTEGER :: ikv(3,20), dikv(3,3), ii
349- REAL(8) :: x, y, z, xv(3)
350- !
351- ! Search nearest neighbor grid points.
348+ INTEGER :: ii
349+ REAL(8) :: xv(3)
350+ integer :: ikv0(3), ikv1(3), i1, i2, i3
352351 !
353352 xv(1:3) = kvec(1:3) * DBLE(ng(1:3))
354- ikv(1:3,1) = NINT(xv(1:3))
355- dikv(1:3,1:3) = 0
356- DO ii = 1, 3
357- dikv(ii,ii) = ikv(ii,1) - FLOOR(xv(ii))
358- dikv(ii,ii) = 1 - 2 * dikv(ii,ii)
359- END DO
360- xv(1:3) = ABS(xv(1:3) - DBLE(ikv(1:3,1)))
361- x = xv(1)
362- y = xv(2)
363- z = xv(3)
353+ ikv0(1:3) = FLOOR(xv(1:3))
354+ xv(1:3) = xv(1:3) - DBLE(ikv0(1:3))
364355 !
365- ikv(1:3, 2) = ikv(1:3,1) + dikv(1:3,1)
366- ikv(1:3, 3) = ikv(1:3,1) + dikv(1:3,2)
367- ikv(1:3, 4) = ikv(1:3,1) + dikv(1:3,3)
368- !
369- IF(nintp == 4) THEN
370- !
371- wintp(1) = 1d0 - x - y - z
372- wintp(2) = x
373- wintp(3) = y
374- wintp(4) = z
375- !
376- ELSE
377- !
378- ikv(1:3, 5) = ikv(1:3,1) + SUM(dikv(1:3,1:3), 2)
379- !
380- ikv(1:3, 6) = ikv(1:3,1) - dikv(1:3,1)
381- ikv(1:3, 7) = ikv(1:3,1) - dikv(1:3,2)
382- ikv(1:3, 8) = ikv(1:3,1) - dikv(1:3,3)
383- !
384- ikv(1:3, 9) = ikv(1:3,1) + 2*dikv(1:3,1)
385- ikv(1:3,10) = ikv(1:3,1) + 2*dikv(1:3,2)
386- ikv(1:3,11) = ikv(1:3,1) + 2*dikv(1:3,3)
387- !
388- ikv(1:3,12) = ikv(1:3,1) + dikv(1:3,2) + dikv(1:3,3)
389- ikv(1:3,13) = ikv(1:3,1) + dikv(1:3,3) + dikv(1:3,1)
390- ikv(1:3,14) = ikv(1:3,1) + dikv(1:3,1) + dikv(1:3,2)
391- !
392- ikv(1:3,15) = ikv(1:3,1) - dikv(1:3,1) + dikv(1:3,3)
393- ikv(1:3,16) = ikv(1:3,1) - dikv(1:3,2) + dikv(1:3,1)
394- ikv(1:3,17) = ikv(1:3,1) - dikv(1:3,3) + dikv(1:3,2)
395- !
396- ikv(1:3,18) = ikv(1:3,1) + dikv(1:3,1) - dikv(1:3,3)
397- ikv(1:3,19) = ikv(1:3,1) + dikv(1:3,2) - dikv(1:3,1)
398- ikv(1:3,20) = ikv(1:3,1) + dikv(1:3,3) - dikv(1:3,2)
399- !
400- wintp( 1) = ( (x - 2d0)*(x - 1d0)*(1d0 + x) &
401- & + (y - 2d0)*(y - 1d0)*(1d0 + y) &
402- & + (z - 2d0)*(z - 1d0)*(1d0 + z) &
403- & + 2d0*(x*y + y*z + z*x)*(x + y + z - 1d0) &
404- & - 8d0*x*y*z - 4d0) * 0.5d0
405- wintp( 2) = x * ( 2d0 + x*(1d0 - x - y - z) &
406- & + y*(1d0 - 2d0*y + z) &
407- & + z*(1d0 - 2d0*z + y)) * 0.5d0
408- wintp( 3) = y * ( 2d0 + y*(1d0 - x - y - z) &
409- & + x*(1d0 - 2d0*x + z) &
410- & + z*(1d0 - 2d0*z + x)) * 0.5d0
411- wintp( 4) = z * ( 2d0 + z*(1d0 - x - y - z) &
412- & + y*(1d0 - 2d0*y + x) &
413- & + x*(1d0 - 2d0*x + y)) * 0.5d0
414- wintp( 5) = x * y * z
415- wintp( 6) = x * (1d0 - x) * ( x + 3d0*y + 3d0*z - 2d0) / 6d0
416- wintp( 7) = y * (1d0 - y) * (3d0*x + y + 3d0*z - 2d0) / 6d0
417- wintp( 8) = z * (1d0 - z) * (3d0*x + 3d0*y + z - 2d0) / 6d0
418- wintp( 9) = x * (x - 1d0) * (x + 1d0) / 6d0
419- wintp(10) = y * (y - 1d0) * (y + 1d0) / 6d0
420- wintp(11) = z * (z - 1d0) * (z + 1d0) / 6d0
421- wintp(12) = y * z * (y + z - 2d0 * x) * 0.5d0
422- wintp(13) = z * x * (z + x - 2d0 * y) * 0.5d0
423- wintp(14) = x * y * (x + y - 2d0 * z) * 0.5d0
424- wintp(15) = x * z * (x - 1d0) * 0.5d0
425- wintp(16) = x * y * (y - 1d0) * 0.5d0
426- wintp(17) = y * z * (z - 1d0) * 0.5d0
427- wintp(18) = x * z * (z - 1d0) * 0.5d0
428- wintp(19) = x * y * (x - 1d0) * 0.5d0
429- wintp(20) = y * z * (y - 1d0) * 0.5d0
430- !
431- END IF
432- !
433- DO ii = 1, nintp
434- ikv(1:3,ii) = MODULO(ikv(1:3,ii), ng(1:3))
435- kintp(ii) = 1 + ikv(1,ii) + ng(1) * ikv(2,ii) + ng(1) * ng(2) * ikv(3,ii)
356+ ii = 0
357+ DO i1 = 0, 1
358+ DO i2 = 0, 1
359+ DO i3 = 0, 1
360+ !
361+ ii = ii + 1
362+ !
363+ ikv1(1:3) = ikv0(1:3) + (/i1, i2, i3/)
364+ ikv1(1:3) = MODULO(ikv1(1:3), ng(1:3))
365+ kintp(ii) = 1 + ikv1(1) + ng(1) * ikv1(2) + ng(1) * ng(2) * ikv1(3)
366+ !
367+ wintp(ii) = xv(1)**i1 * (1.0d0 - xv(1))**(1 - i1) &
368+ & * xv(2)**i2 * (1.0d0 - xv(2))**(1 - i2) &
369+ & * xv(3)**i3 * (1.0d0 - xv(3))**(1 - i3)
370+ !
371+ END DO
372+ END DO
436373 END DO
437374 !
438375 END SUBROUTINE libtetrabz_interpol_indx
--- a/src/libtetrabz_dbldelta_mod.F90
+++ b/src/libtetrabz_dbldelta_mod.F90
@@ -43,13 +43,11 @@ SUBROUTINE libtetrabz_dbldelta(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm) BIND(
4343 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
4444 !
4545 LOGICAL :: linterpol
46- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
46+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
4747 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
48- REAL(8) :: wlsm(4,20), wintp(1,20)
48+ REAL(8) :: wlsm(4,20), wintp(1,8)
4949 REAL(8),ALLOCATABLE :: wghtd(:,:,:), kvec(:,:)
5050 !
51- nintp = 16 * ltetra - 12
52- !
5351 IF(PRESENT(comm)) THEN
5452 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
5553 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -67,9 +65,9 @@ SUBROUTINE libtetrabz_dbldelta(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm) BIND(
6765 !
6866 wght(1:nb*nb,1:PRODUCT(ngw(1:3))) = 0d0
6967 DO ik = 1, nk_local
70- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
71- wght(1:nb*nb,kintp(1:nintp)) = wght(1:nb*nb, kintp(1:nintp)) &
72- & + MATMUL(wghtd(1:nb*nb,1:1,ik), wintp(1:1,1:nintp))
68+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
69+ wght(1:nb*nb,kintp(1:8)) = wght(1:nb*nb, kintp(1:8)) &
70+ & + MATMUL(wghtd(1:nb*nb,1:1,ik), wintp(1:1,1:8))
7371 END DO ! ik = 1, nk_local
7472 DEALLOCATE(wghtd, kvec)
7573 !
--- a/src/libtetrabz_dblstep_mod.F90
+++ b/src/libtetrabz_dblstep_mod.F90
@@ -43,13 +43,11 @@ SUBROUTINE libtetrabz_dblstep(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm) BIND(C
4343 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
4444 !
4545 LOGICAL :: linterpol
46- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
46+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
4747 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
48- REAL(8) :: wlsm(4,20), wintp(1,20)
48+ REAL(8) :: wlsm(4,20), wintp(1,8)
4949 REAL(8),ALLOCATABLE :: wghtd(:,:,:), kvec(:,:)
5050 !
51- nintp = 16 * ltetra - 12
52- !
5351 IF(PRESENT(comm)) THEN
5452 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
5553 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -67,9 +65,9 @@ SUBROUTINE libtetrabz_dblstep(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm) BIND(C
6765 !
6866 wght(1:nb*nb,1:PRODUCT(ngw(1:3))) = 0d0
6967 DO ik = 1, nk_local
70- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
71- wght(1:nb*nb,kintp(1:nintp)) = wght(1:nb*nb, kintp(1:nintp)) &
72- & + MATMUL(wghtd(1:nb*nb,1:1,ik), wintp(1:1,1:nintp))
68+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
69+ wght(1:nb*nb,kintp(1:8)) = wght(1:nb*nb, kintp(1:8)) &
70+ & + MATMUL(wghtd(1:nb*nb,1:1,ik), wintp(1:1,1:8))
7371 END DO ! ik = 1, nk_local
7472 DEALLOCATE(wghtd, kvec)
7573 !
--- a/src/libtetrabz_dos_mod.F90
+++ b/src/libtetrabz_dos_mod.F90
@@ -43,13 +43,11 @@ SUBROUTINE libtetrabz_dos(ltetra,bvec,nb,nge,eig,ngw,wght,ne,e0,comm) BIND(C)
4343 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
4444 !
4545 LOGICAL :: linterpol
46- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
46+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
4747 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
48- REAL(8) :: wlsm(4,20), wintp(1,20)
48+ REAL(8) :: wlsm(4,20), wintp(1,8)
4949 REAL(8),ALLOCATABLE :: wghtd(:,:,:), kvec(:,:)
5050 !
51- nintp = 16 * ltetra - 12
52- !
5351 IF(PRESENT(comm)) THEN
5452 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
5553 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -67,9 +65,9 @@ SUBROUTINE libtetrabz_dos(ltetra,bvec,nb,nge,eig,ngw,wght,ne,e0,comm) BIND(C)
6765 !
6866 wght(1:ne*nb,1:PRODUCT(ngw(1:3))) = 0d0
6967 DO ik = 1, nk_local
70- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
71- wght(1:ne*nb,kintp(1:nintp)) = wght(1:ne*nb, kintp(1:nintp)) &
72- & + MATMUL(wghtd(1:ne*nb,1:1,ik), wintp(1:1,1:nintp))
68+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
69+ wght(1:ne*nb,kintp(1:8)) = wght(1:ne*nb, kintp(1:8)) &
70+ & + MATMUL(wghtd(1:ne*nb,1:1,ik), wintp(1:1,1:8))
7371 END DO ! ik = 1, nk_local
7472 DEALLOCATE(wghtd, kvec)
7573 !
@@ -97,13 +95,11 @@ SUBROUTINE libtetrabz_intdos(ltetra,bvec,nb,nge,eig,ngw,wght,ne,e0,comm) BIND(C)
9795 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
9896 !
9997 LOGICAL :: linterpol
100- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
98+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
10199 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
102- REAL(8) :: wlsm(4,20), wintp(1,20)
100+ REAL(8) :: wlsm(4,20), wintp(1,8)
103101 REAL(8),ALLOCATABLE :: wghtd(:,:,:), kvec(:,:)
104102 !
105- nintp = 16 * ltetra - 12
106- !
107103 IF(PRESENT(comm)) THEN
108104 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
109105 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -121,9 +117,9 @@ SUBROUTINE libtetrabz_intdos(ltetra,bvec,nb,nge,eig,ngw,wght,ne,e0,comm) BIND(C)
121117 !
122118 wght(1:ne*nb,1:PRODUCT(ngw(1:3))) = 0d0
123119 DO ik = 1, nk_local
124- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
125- wght(1:ne*nb,kintp(1:nintp)) = wght(1:ne*nb, kintp(1:nintp)) &
126- & + MATMUL(wghtd(1:ne*nb,1:1,ik), wintp(1:1,1:nintp))
120+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
121+ wght(1:ne*nb,kintp(1:8)) = wght(1:ne*nb, kintp(1:8)) &
122+ & + MATMUL(wghtd(1:ne*nb,1:1,ik), wintp(1:1,1:8))
127123 END DO ! ik = 1, nk_local
128124 DEALLOCATE(wghtd, kvec)
129125 !
--- a/src/libtetrabz_fermigr_mod.F90
+++ b/src/libtetrabz_fermigr_mod.F90
@@ -43,13 +43,11 @@ SUBROUTINE libtetrabz_fermigr(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,ne,e0,comm)
4343 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
4444 !
4545 LOGICAL :: linterpol
46- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
46+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
4747 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
48- REAL(8) :: wlsm(4,20), wintp(1,20)
48+ REAL(8) :: wlsm(4,20), wintp(1,8)
4949 REAL(8),ALLOCATABLE :: wghtd(:,:,:), kvec(:,:)
5050 !
51- nintp = 16 * ltetra - 12
52- !
5351 IF(PRESENT(comm)) THEN
5452 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
5553 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -67,9 +65,9 @@ SUBROUTINE libtetrabz_fermigr(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,ne,e0,comm)
6765 !
6866 wght(1:ne*nb*nb,1:PRODUCT(ngw(1:3))) = 0d0
6967 DO ik = 1, nk_local
70- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
71- wght(1:ne*nb*nb,kintp(1:nintp)) = wght(1:ne*nb*nb, kintp(1:nintp)) &
72- & + MATMUL(wghtd(1:ne*nb*nb,1:1,ik), wintp(1:1,1:nintp))
68+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
69+ wght(1:ne*nb*nb,kintp(1:8)) = wght(1:ne*nb*nb, kintp(1:8)) &
70+ & + MATMUL(wghtd(1:ne*nb*nb,1:1,ik), wintp(1:1,1:8))
7371 END DO ! ik = 1, nk_local
7472 DEALLOCATE(wghtd, kvec)
7573 !
--- a/src/libtetrabz_occ_mod.F90
+++ b/src/libtetrabz_occ_mod.F90
@@ -43,13 +43,11 @@ SUBROUTINE libtetrabz_occ(ltetra,bvec,nb,nge,eig,ngw,wght,comm) BIND(C)
4343 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
4444 !
4545 LOGICAL :: linterpol
46- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
46+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
4747 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
48- REAL(8) :: wlsm(4,20), wintp(1,20)
48+ REAL(8) :: wlsm(4,20), wintp(1,8)
4949 REAL(8),ALLOCATABLE :: wghtd(:,:,:), kvec(:,:)
5050 !
51- nintp = 16 * ltetra - 12
52- !
5351 IF(PRESENT(comm)) THEN
5452 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
5553 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -67,9 +65,9 @@ SUBROUTINE libtetrabz_occ(ltetra,bvec,nb,nge,eig,ngw,wght,comm) BIND(C)
6765 !
6866 wght(1:nb,1:PRODUCT(ngw(1:3))) = 0d0
6967 DO ik = 1, nk_local
70- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
71- wght(1:nb,kintp(1:nintp)) = wght(1:nb, kintp(1:nintp)) &
72- & + MATMUL(wghtd(1:nb,1:1,ik), wintp(1:1,1:nintp))
68+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
69+ wght(1:nb,kintp(1:8)) = wght(1:nb, kintp(1:8)) &
70+ & + MATMUL(wghtd(1:nb,1:1,ik), wintp(1:1,1:8))
7371 END DO ! ik = 1, nk_local
7472 DEALLOCATE(wghtd, kvec)
7573 !
@@ -99,16 +97,14 @@ SUBROUTINE libtetrabz_fermieng(ltetra,bvec,nb,nge,eig,ngw,wght,ef,nelec,comm) BI
9997 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
10098 !
10199 LOGICAL :: linterpol
102- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
100+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
103101 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
104- REAL(8) :: wlsm(4,20), wintp(1,20)
102+ REAL(8) :: wlsm(4,20), wintp(1,8)
105103 REAL(8),ALLOCATABLE :: wghtd(:,:,:), kvec(:,:)
106104 !
107105 INTEGER :: iter, maxiter = 300
108106 REAL(8) :: elw, eup, eps= 1d-10, sumkmid
109107 !
110- nintp = 16 * ltetra - 12
111- !
112108 IF(PRESENT(comm)) THEN
113109 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
114110 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -161,9 +157,9 @@ SUBROUTINE libtetrabz_fermieng(ltetra,bvec,nb,nge,eig,ngw,wght,ef,nelec,comm) BI
161157 !
162158 wght(1:nb,1:PRODUCT(ngw(1:3))) = 0d0
163159 DO ik = 1, nk_local
164- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
165- wght(1:nb,kintp(1:nintp)) = wght(1:nb, kintp(1:nintp)) &
166- & + MATMUL(wghtd(1:nb,1:1,ik), wintp(1:1,1:nintp))
160+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
161+ wght(1:nb,kintp(1:8)) = wght(1:nb, kintp(1:8)) &
162+ & + MATMUL(wghtd(1:nb,1:1,ik), wintp(1:1,1:8))
167163 END DO ! ik = 1, nk_local
168164 DEALLOCATE(wghtd, kvec)
169165 !
--- a/src/libtetrabz_polcmplx_mod.F90
+++ b/src/libtetrabz_polcmplx_mod.F90
@@ -44,14 +44,12 @@ SUBROUTINE libtetrabz_polcmplx(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,ne,e0,comm)
4444 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
4545 !
4646 LOGICAL :: linterpol
47- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
47+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
4848 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
49- REAL(8) :: wlsm(4,20), wintp(1,20)
49+ REAL(8) :: wlsm(4,20), wintp(1,8)
5050 REAL(8),ALLOCATABLE :: kvec(:,:)
5151 COMPLEX(8),ALLOCATABLE :: wghtd(:,:,:)
5252 !
53- nintp = 16 * ltetra - 12
54- !
5553 IF(PRESENT(comm)) THEN
5654 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
5755 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -69,9 +67,9 @@ SUBROUTINE libtetrabz_polcmplx(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,ne,e0,comm)
6967 !
7068 wght(1:ne*nb*nb,1:PRODUCT(ngw(1:3))) = 0d0
7169 DO ik = 1, nk_local
72- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
73- wght(1:ne*nb*nb,kintp(1:nintp)) = wght(1:ne*nb*nb, kintp(1:nintp)) &
74- & + MATMUL(wghtd(1:ne*nb*nb,1:1,ik), wintp(1:1,1:nintp))
70+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
71+ wght(1:ne*nb*nb,kintp(1:8)) = wght(1:ne*nb*nb, kintp(1:8)) &
72+ & + MATMUL(wghtd(1:ne*nb*nb,1:1,ik), wintp(1:1,1:8))
7573 END DO ! ik = 1, nk_local
7674 DEALLOCATE(wghtd, kvec)
7775 !
--- a/src/libtetrabz_polstat_mod.F90
+++ b/src/libtetrabz_polstat_mod.F90
@@ -43,13 +43,11 @@ SUBROUTINE libtetrabz_polstat(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm) BIND(C
4343 INTEGER(C_INT),INTENT(IN),OPTIONAL :: comm
4444 !
4545 LOGICAL :: linterpol
46- INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(20), nintp
46+ INTEGER :: nt_local, nk_local, nkBZ, ik, kintp(8)
4747 INTEGER,ALLOCATABLE :: ik_global(:,:), ik_local(:,:)
48- REAL(8) :: wlsm(4,20), wintp(1,20)
48+ REAL(8) :: wlsm(4,20), wintp(1,8)
4949 REAL(8),ALLOCATABLE :: wghtd(:,:,:), kvec(:,:)
5050 !
51- nintp = 16 * ltetra - 12
52- !
5351 IF(PRESENT(comm)) THEN
5452 CALL libtetrabz_initialize(ltetra,nge,ngw,bvec,linterpol,wlsm,nk_local,&
5553 & nt_local,nkBZ,ik_global,ik_local,kvec,comm)
@@ -67,9 +65,9 @@ SUBROUTINE libtetrabz_polstat(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm) BIND(C
6765 !
6866 wght(1:nb*nb,1:PRODUCT(ngw(1:3))) = 0d0
6967 DO ik = 1, nk_local
70- CALL libtetrabz_interpol_indx(nintp,ngw,kvec(1:3,ik),kintp,wintp)
71- wght(1:nb*nb,kintp(1:nintp)) = wght(1:nb*nb, kintp(1:nintp)) &
72- & + MATMUL(wghtd(1:nb*nb,1:1,ik), wintp(1:1,1:nintp))
68+ CALL libtetrabz_interpol_indx(ngw,kvec(1:3,ik),kintp,wintp)
69+ wght(1:nb*nb,kintp(1:8)) = wght(1:nb*nb, kintp(1:8)) &
70+ & + MATMUL(wghtd(1:nb*nb,1:1,ik), wintp(1:1,1:8))
7371 END DO ! ik = 1, nk_local
7472 DEALLOCATE(wghtd, kvec)
7573 !
--- a/test/test2_16_8.sh
+++ b/test/test2_16_8.sh
@@ -9,95 +9,95 @@ ${MPIRUN} ./test.x < test.in > test.out
99 cat > test.out.ref <<EOF
1010 # Ideal Result
1111 # libtetrabz_occ
12- 0.25133E+01 0.25139E+01
13- 0.44429E+00 0.44363E+00
12+ 0.25133E+01 0.27342E+01
13+ 0.44429E+00 0.52281E+00
1414
1515 # libtetrabz_fermieng
1616 0.50000E+00 0.50001E+00
17- 0.25133E+01 0.25141E+01
18- 0.44429E+00 0.44369E+00
17+ 0.25133E+01 0.27344E+01
18+ 0.44429E+00 0.52287E+00
1919
2020 # libtetrabz_dos
21- 0.10053E+00 0.10734E+00
22- 0.80425E+00 0.82145E+00
23- 0.27143E+01 0.27056E+01
24- 0.64340E+01 0.64356E+01
25- 0.12566E+02 0.12523E+02
21+ 0.10053E+00 0.25701E+00
22+ 0.80425E+00 0.10796E+01
23+ 0.27143E+01 0.31247E+01
24+ 0.64340E+01 0.69397E+01
25+ 0.12566E+02 0.13222E+02
2626 0.00000E+00 0.00000E+00
2727 0.00000E+00 0.00000E+00
2828 0.00000E+00 0.00000E+00
29- 0.65827E+00 0.68276E+00
30- 0.44429E+01 0.44352E+01
29+ 0.65827E+00 0.93716E+00
30+ 0.44429E+01 0.49256E+01
3131
3232 # libtetrabz_intdos
33- 0.80425E-03 0.82841E-03
34- 0.25736E-01 0.25976E-01
35- 0.19543E+00 0.19585E+00
36- 0.82355E+00 0.82435E+00
37- 0.25133E+01 0.25139E+01
33+ 0.80425E-03 0.22258E-02
34+ 0.25736E-01 0.40841E-01
35+ 0.19543E+00 0.24277E+00
36+ 0.82355E+00 0.93854E+00
37+ 0.25133E+01 0.27342E+01
3838 0.00000E+00 0.00000E+00
3939 0.00000E+00 0.00000E+00
4040 0.00000E+00 0.00000E+00
41- 0.18431E-01 0.18434E-01
42- 0.44429E+00 0.44363E+00
41+ 0.18431E-01 0.30728E-01
42+ 0.44429E+00 0.52281E+00
4343
4444 # libtetrabz_dblstep
45- 0.48106E+00 0.48116E+00
45+ 0.48106E+00 0.51756E+00
4646 0.00000E+00 0.00000E+00
47- 0.12428E+00 0.12419E+00
47+ 0.12428E+00 0.14246E+00
4848 0.00000E+00 0.00000E+00
4949
5050 # libtetrabz_dbldelta
51- 0.62832E+01 0.62758E+01
51+ 0.62832E+01 0.66298E+01
5252 0.00000E+00 0.00000E+00
53- 0.31416E+01 0.31276E+01
53+ 0.31416E+01 0.35252E+01
5454 0.00000E+00 0.00000E+00
5555
5656 # libtetrabz_polstat
57- 0.38431E+01 0.38434E+01
58- 0.50265E+01 0.50278E+01
59- 0.10489E+01 0.10445E+01
60- 0.17772E+01 0.17745E+01
57+ 0.38431E+01 0.41364E+01
58+ 0.50265E+01 0.54685E+01
59+ 0.10489E+01 0.12224E+01
60+ 0.17772E+01 0.20912E+01
6161
6262 # libtetrabz_fermigr
63- 0.13963E+01 0.13962E+01
64- 0.15696E+01 0.15704E+01
65- 0.14726E+01 0.14696E+01
63+ 0.13963E+01 0.15389E+01
64+ 0.15696E+01 0.17770E+01
65+ 0.14726E+01 0.16065E+01
6666 0.00000E+00 0.00000E+00
6767 0.00000E+00 0.00000E+00
6868 0.00000E+00 0.00000E+00
69- 0.39262E+00 0.39179E+00
70- 0.34535E+00 0.34688E+00
69+ 0.39262E+00 0.47167E+00
70+ 0.34535E+00 0.38944E+00
7171 0.00000E+00 0.00000E+00
7272 0.00000E+00 0.00000E+00
7373 0.00000E+00 0.00000E+00
7474 0.00000E+00 0.00000E+00
7575
7676 # libtetrabz_polcmplx
77- -0.83824E+00 -0.83659E+00
78- -0.73420E+00 -0.73974E+00
79- 0.27039E+00 0.27046E+00
80- -0.77191E+00 -0.77057E+00
81- 0.97100E+00 0.97111E+00
82- 0.30279E+00 0.30382E+00
83- -0.11600E+01 -0.11603E+01
84- -0.77332E+00 -0.77351E+00
85- 0.29568E+00 0.29575E+00
86- -0.11827E+01 -0.11830E+01
87- 0.15080E+01 0.15083E+01
88- 0.50265E+00 0.50278E+00
89- -0.13077E+00 -0.13063E+00
90- -0.87431E-01 -0.88020E-01
91- 0.30122E-01 0.30298E-01
92- -0.13535E+00 -0.13475E+00
93- 0.17888E+00 0.17820E+00
94- 0.64232E-01 0.64093E-01
95- -0.19139E+00 -0.19110E+00
96- -0.10936E+00 -0.10920E+00
97- 0.27341E-01 0.27301E-01
98- -0.21873E+00 -0.21841E+00
99- 0.30641E+00 0.30596E+00
100- 0.12256E+00 0.12238E+00
77+ -0.83824E+00 -0.90780E+00
78+ -0.73420E+00 -0.80092E+00
79+ 0.27039E+00 0.29335E+00
80+ -0.77191E+00 -0.83629E+00
81+ 0.97100E+00 0.10534E+01
82+ 0.30279E+00 0.32932E+00
83+ -0.11600E+01 -0.12620E+01
84+ -0.77332E+00 -0.84130E+00
85+ 0.29568E+00 0.32168E+00
86+ -0.11827E+01 -0.12867E+01
87+ 0.15080E+01 0.16405E+01
88+ 0.50265E+00 0.54685E+00
89+ -0.13077E+00 -0.15323E+00
90+ -0.87431E-01 -0.10300E+00
91+ 0.30122E-01 0.35394E-01
92+ -0.13535E+00 -0.15828E+00
93+ 0.17888E+00 0.20943E+00
94+ 0.64232E-01 0.75401E-01
95+ -0.19139E+00 -0.22521E+00
96+ -0.10936E+00 -0.12869E+00
97+ 0.27341E-01 0.32173E-01
98+ -0.21873E+00 -0.25738E+00
99+ 0.30641E+00 0.36056E+00
100+ 0.12256E+00 0.14422E+00
101101
102102 EOF
103103
旧リポジトリブラウザで表示