TABLE OF CONTENTS


ABINIT/m_GreenHyboffdiag [ Modules ]

[ Top ] [ Modules ]

NAME

  m_GreenHyboffdiag

FUNCTION

  Manage a green function for one orbital

COPYRIGHT

  Copyright (C) 2013 ABINIT group (B.Amadon, J. Denier and J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

28 #include "defs.h"
29 MODULE m_GreenHyboffdiag
30 
31  USE m_global
32  USE m_MatrixHyb
33  USE m_Vector
34  USE m_VectorInt
35  USE m_ListCdagC
36  USE m_MapHyb
37 #ifdef HAVE_MPI2
38  USE mpi
39 #endif
40  
41  IMPLICIT NONE
42 
43 
44  public ::  GreenHyboffdiag_init
45  public ::  GreenHyboffdiag_reset
46  public ::  GreenHyboffdiag_clear
47  public ::  GreenHyboffdiag_setOperW
48  public ::  GreenHyboffdiag_measHybrid
49  public ::  GreenHyboffdiag_getHybrid
50  public ::  GreenHyboffdiag_setN
51  public ::  GreenHyboffdiag_setMuD1
52  public ::  GreenHyboffdiag_setMoments
53  public ::  GreenHyboffdiag_backFourier
54  public ::  GreenHyboffdiag_forFourier
55  public ::  GreenHyboffdiag_print
56  public ::  GreenHyboffdiag_destroy
57  public ::  nfourier3

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_backFourier [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_backFourier

FUNCTION

  perform back fourier transform

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green
  dvgc=divergence parameter

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

1168 SUBROUTINE GreenHyboffdiag_backFourier(op,dvgc,func,hybri_limit,opt_hybri_limit)
1169 
1170  use m_fstrings,       only : int2char4
1171 
1172 !This section has been created automatically by the script Abilint (TD).
1173 !Do not modify the following lines by hand.
1174 #undef ABI_FUNC
1175 #define ABI_FUNC 'GreenHyboffdiag_backFourier'
1176 !End of the abilint section
1177 
1178 
1179 #ifdef HAVE_MPI1
1180 include 'mpif.h'
1181 #endif
1182 !Arguments ------------------------------------
1183   TYPE(GreenHyboffdiag)            , INTENT(INOUT) :: op
1184   DOUBLE PRECISION, OPTIONAL, INTENT(IN   ) :: dvgc
1185   CHARACTER(len=5)  ,OPTIONAL, INTENT(IN) :: func
1186   COMPLEX(KIND=8), DIMENSION(op%nflavors,op%nflavors), OPTIONAL, INTENT(IN) :: hybri_limit
1187   INTEGER, OPTIONAL, INTENT(IN) :: opt_hybri_limit
1188 !Local variables ------------------------------
1189   INTEGER :: itau
1190   INTEGER :: iomega
1191   INTEGER :: omegaSamples
1192   INTEGER :: tauSamples
1193   INTEGER :: tauBegin
1194   INTEGER :: tauEnd
1195   INTEGER :: delta
1196   INTEGER :: residu
1197   INTEGER :: iflavor1
1198   INTEGER :: iflavor2,unitnb,unitnb1
1199   INTEGER, ALLOCATABLE, DIMENSION(:) :: counts
1200   INTEGER, ALLOCATABLE, DIMENSION(:) :: displs
1201   DOUBLE PRECISION :: A,AA ! Correction factor
1202   COMPLEX(KIND=8) :: B,BB ! Correction factor
1203   COMPLEX(KIND=8) :: C,CC ! Correction factor
1204   DOUBLE PRECISION :: inv_beta
1205   DOUBLE PRECISION :: pi_invBeta
1206   DOUBLE PRECISION :: two_invBeta
1207   DOUBLE PRECISION :: minusDt
1208   DOUBLE PRECISION :: minusOmegaTau
1209   DOUBLE PRECISION :: omegaa
1210   DOUBLE PRECISION :: minusTau
1211   DOUBLE PRECISION :: sumTerm
1212   DOUBLE PRECISION :: pi
1213   DOUBLE PRECISION :: twoPi
1214   DOUBLE PRECISION :: correction
1215   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: Domega
1216   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: A_omega
1217   COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:) :: C_omega
1218   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: opertau
1219   CHARACTER(len=5) :: funct
1220   character(len=4) :: tag_proc
1221   character(len=30) :: tmpfil
1222 
1223   IF ( op%set .EQV. .FALSE. ) &
1224     CALL ERROR("GreenHyboffdiag_backFourier : Uninitialized GreenHyboffdiag structure")
1225   IF ( op%setW .EQV. .FALSE. ) &
1226     CALL ERROR("GreenHyboffdiag_backFourier : no G(iw)")
1227   
1228   funct="hybri"
1229   if(present(func)) funct=func
1230 !sui!write(6,*) funct
1231   inv_beta     = op%inv_beta
1232   two_invBeta  = 2.d0 * inv_beta
1233   minusDt      = - op%delta_t
1234   omegaSamples = op%Wmax
1235   tauSamples   = op%samples-1
1236   pi         = ACOS(-1.d0)
1237   twoPi        = 2.d0 * pi
1238   pi_invBeta = pi * inv_beta
1239 !sui!write(6,*) "omegaSamples",omegaSamples
1240   MALLOC(Domega,(1:omegaSamples))
1241   MALLOC(A_omega,(1:omegaSamples))
1242   MALLOC(C_omega,(1:omegaSamples))
1243   IF ( op%rank .EQ. 0 ) THEN
1244     !DO iflavor1 = 1, op%nflavors
1245     !  DO iflavor2 = 1, op%nflavors
1246     !    write(22236,*) "#",iflavor1,iflavor2
1247     !    do  iomega=1,op%Wmax
1248     !      write(22236,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2))
1249     !    enddo
1250     !    write(22236,*) 
1251     !  ENDDO
1252     !ENDDO
1253   ENDIF
1254 
1255   op%oper = 0.d0
1256 
1257   DO iflavor1 = 1, op%nflavors
1258     DO iflavor2 = 1, op%nflavors
1259   ! --  compute limit of function G*(i\omega_n)
1260       if(funct=="hybri") then
1261         IF ( PRESENT(dvgc) ) THEN
1262           A = dvgc
1263         ELSE
1264           A = AIMAG(op%oper_w(omegaSamples,iflavor1,iflavor2))&! A = \lim_\infty Imag(G(iwn))*(wn)
1265             *(2.d0*DBLE(omegaSamples)-1.d0) * pi_invBeta
1266           AA = AIMAG(op%oper_w(omegaSamples-10,iflavor1,iflavor2))&! A = \lim_\infty Imag(G(iwn))*(wn)
1267             *(2.d0*DBLE(omegaSamples-10)-1.d0) * pi_invBeta
1268           B = op%oper_w(omegaSamples,iflavor1,iflavor2)&! A = \lim_\infty (G(iwn))*(iwn)
1269             *(2.d0*DBLE(omegaSamples)-1.d0) *pi_invBeta*cmplx(0.d0,1.d0,kind=8)
1270           BB = op%oper_w(omegaSamples-10,iflavor1,iflavor2)&! A = \lim_\infty (G(iwn))*(iwn)
1271             *(2.d0*DBLE(omegaSamples-10)-1.d0) *pi_invBeta*cmplx(0.d0,1.d0,kind=8)
1272         !sui!write(6,*) "B=",iflavor1,iflavor2,B,BB
1273         END IF
1274       else if(iflavor1==iflavor2.and.funct=="green") then
1275         A = -1.d0
1276       else if(iflavor1/=iflavor2.and.funct=="green") then
1277         A = 0.d0
1278       endif
1279     !sui!write(6,*) "A=",iflavor1,iflavor2,A,AA
1280       C=cmplx(-A,0.d0,kind=8)
1281       if(present(hybri_limit)) then
1282         if(present(opt_hybri_limit)) then
1283           if(opt_hybri_limit==1) C= (hybri_limit(iflavor1,iflavor2))
1284         !sui!write(6,*) "C=                         ",C
1285         endif
1286       endif 
1287 
1288 
1289 
1290   ! --  correction on G(tau=0) is thus
1291       correction = -C*0.5d0
1292     
1293   ! --  built frequency mesh
1294       Domega = (/ ((2.d0 * DBLE(iomega) - 1.d0)*pi_invbeta, iomega=1, omegaSamples) /)
1295 
1296   ! --  built asymptotic function
1297       !if(present(hybri_limit)) then
1298       C_omega = C / (Domega*cmplx(0.d0,1.d0,kind=8))
1299       !else
1300       !  A_omega = A / Domega
1301       !  C_omega=cmplx(0.d0,A_omega,kind=8)
1302       !endif
1303       !write(6,*) "AC 1",A_omega(2),C_omega(2)
1304 
1305       IF ( op%rank .EQ. 0 ) THEN
1306        ! write(236,*) "#",iflavor1,iflavor2
1307        ! write(237,*) "#",iflavor1,iflavor2
1308        ! write(2236,*) "#",iflavor1,iflavor2
1309        ! write(2237,*) "#",iflavor1,iflavor2
1310        ! write(238,*) "#",iflavor1,iflavor2
1311        ! do  iomega=1,op%Wmax
1312        !   write(2236,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2))
1313        !   write(2237,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,imag(op%oper_w(iomega,iflavor1,iflavor2))
1314        !   write(236,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2)-C_omega(iomega))
1315        !   write(237,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,imag(op%oper_w(iomega,iflavor1,iflavor2)-C_omega(iomega))
1316        !   write(238,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,C_omega(iomega)
1317        ! enddo
1318        ! write(236,*) 
1319        ! write(237,*) 
1320        ! write(2236,*) 
1321        ! write(2237,*) 
1322        ! write(238,*) 
1323       ENDIF
1324       IF ( op%rank .EQ. 1 ) THEN
1325        ! write(22360,*) "#",iflavor1,iflavor2
1326        ! do  iomega=1,op%Wmax
1327        !   write(22360,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2))
1328        ! enddo
1329        ! write(22360,*) 
1330       ENDIF
1331       !open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted')
1332       !write(unitnb,*) "#",iflavor1,iflavor2
1333       !do  iomega=1,op%Wmax
1334       !  write(unitnb,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2))
1335       !enddo
1336       !write(unitnb,*) 
1337   ! --  built time mesh
1338       IF (op%have_MPI .EQV. .TRUE.) THEN
1339         delta = tauSamples / op%size
1340         residu = tauSamples - op%size*delta
1341         IF ( op%rank .LT. op%size - residu ) THEN
1342           tauBegin = 1 + op%rank*delta
1343           tauEnd   = (op%rank + 1)*delta
1344         ELSE
1345 !          tauBegin = (op%size-residu)*delta + 1 + (op%rank-op%size+residu)*(delta+1)
1346           tauBegin = 1 + op%rank*(delta + 1) -op%size + residu
1347           tauEnd = tauBegin + delta
1348         END IF
1349         MALLOC(counts,(1:op%size))
1350         MALLOC(displs,(1:op%size))
1351         counts = (/ (delta, iTau=1, op%size-residu), &
1352                     (delta+1, iTau=op%size-residu+1, op%size) /)
1353         displs(1)=0
1354         DO iTau = 2, op%size
1355           displs(iTau) = displs(iTau-1) + counts (iTau-1)
1356         END DO
1357       ELSE
1358         tauBegin = 1
1359         tauEnd   = tauSamples
1360       END IF
1361       MALLOC(opertau,(1:tauSamples+1))
1362       do iomega=1,omegaSamples
1363        ! write(6,*) iomega, imag(op%oper_w(iomega,iflavor1,iflavor2)), A_omega(iomega) ,"#diff"
1364       enddo
1365       unitnb=70000+op%rank
1366       call int2char4(op%rank,tag_proc)
1367       tmpfil = 'counts'//tag_proc
1368      ! open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted')
1369      ! write(unitnb,*) "#",iflavor1,iflavor2
1370      ! do  itau=1,op%size
1371      ! write(unitnb,*)  itau,counts(itau),displs(itau)
1372      ! enddo
1373      ! write(unitnb,*) 
1374 
1375       unitnb=10000+op%rank
1376       call int2char4(op%rank,tag_proc)
1377       tmpfil = 'oper_w'//tag_proc
1378      ! open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted')
1379      ! write(unitnb,*) "#",iflavor1,iflavor2,C
1380      ! ! C_omega et oper_w differents Domega identique. Est ce du a des
1381      !! ! diago differentes   pour chaque procs dans qmc_prep_ctqmc
1382      !! do  iomega=1,op%Wmax
1383      !! write(unitnb,*)  (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2)),C_omega(iomega),Domega(iomega)
1384      ! enddo
1385      ! write(unitnb,*) 
1386 
1387      ! unitnb=40000+op%rank
1388      ! unitnb1=50000+op%rank
1389      ! call int2char4(op%rank,tag_proc)
1390      ! tmpfil = 'tauend'//tag_proc
1391      ! open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted')
1392      ! tmpfil = 'taubegin'//tag_proc
1393      ! open (unit=unitnb1,file=trim(tmpfil),status='unknown',form='formatted')
1394      ! write(unitnb,*) "#",iflavor1,iflavor2
1395      ! write(unitnb1,*) "#",iflavor1,iflavor2
1396 
1397   ! -- compute Fourier transformation
1398       opertau=0.d0
1399       DO itau = tauBegin, tauEnd
1400       !DO itau = max(tauBegin-1,1), tauEnd
1401         minusTau = DBLE(itau -1) * minusDt
1402         DO iomega = 1, omegaSamples
1403           omegaa         = Domega(iomega)
1404           minusOmegaTau = MOD(omegaa*minusTau, TwoPi)
1405           sumTerm       = REAL(( op%oper_w(iomega,iflavor1,iflavor2) &
1406                           -  C_omega(iomega) ) &
1407                           !- CMPLX(0.d0, A_omega(iomega),8) ) &
1408                           * EXP( CMPLX(0.d0, minusOmegaTau, 8)))
1409           opertau(itau)  = opertau(itau) + sumTerm
1410 !         Domega et minusomegatau identique MAIS oper_w different
1411             !write(unitnb,*) iomega,Domega(iomega),real(C_omega(iomega)),imag(C_omega(iomega))
1412           !if(itau==tauend) then
1413           !  write(unitnb,*) iomega, sumTerm,opertau(itau),minusOmegaTau,op%oper_w(iomega,iflavor1,iflavor2),Domega(iomega)
1414           !endif 
1415           !if(itau==max(tauBegin-1,1)) then
1416           !  write(unitnb1,*) iomega, sumTerm,opertau(itau),minusOmegaTau,op%oper_w(iomega,iflavor1,iflavor2),Domega(iomega)
1417           !endif 
1418 
1419         END DO
1420          ! if(itau==tauEnd) write(unitnb,*) 
1421          ! if(itau==max(tauBegin-1,1)) write(unitnb1,*) 
1422 !        if(iflavor1==iflavor2) then
1423           opertau(itau) = correction + two_invBeta*opertau(itau)
1424          ! if(itau==tauend) then
1425          !   write(unitnb,*) "final", opertau(itau),correction
1426          ! endif 
1427          ! if(itau==max(tauBegin-1,1)) then
1428          !   write(unitnb1,*) "final",opertau(itau),correction
1429          ! endif 
1430              !write(66666,*) itau, opertau(itau),correction
1431 !        else
1432 !          opertau(itau) =              &
1433 !             two_invBeta*opertau(itau)
1434 !        endif
1435       END DO
1436       !unitnb=20000+op%rank
1437       !call int2char4(op%rank,tag_proc)
1438       !tmpfil = 'opertau'//tag_proc
1439       !open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted')
1440       !write(unitnb,*) "#",iflavor1,iflavor2,tauBegin,tauEnd
1441       !!do  itau=tauBegin, tauEnd
1442       !do  itau=1,tauSamples
1443       !  write(unitnb,*)    itau,opertau(itau)
1444       !enddo
1445       !write(unitnb,*) 
1446       !opertau(tauBegin-1)=0.d0
1447       !opertau(tauEnd+1)=0.d0
1448              !write(66666,*)
1449   
1450   ! -- Gather
1451       IF ( op%have_MPI .EQV. .TRUE. ) THEN
1452 ! rassembler les resultats
1453 #ifdef HAVE_MPI
1454         CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DATATYPE_NULL, &
1455                           opertau, counts, displs, &
1456                           MPI_DOUBLE_PRECISION, op%MY_COMM, residu)
1457 #endif
1458         FREE(counts)
1459         FREE(displs)
1460       END IF
1461      ! unitnb=30000+op%rank
1462      ! call int2char4(op%rank,tag_proc)
1463      ! tmpfil = 'opertau_MPI_'//tag_proc
1464      ! open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted')
1465      ! write(unitnb,*) "#",iflavor1,iflavor2
1466      ! do  itau=tauBegin, tauEnd
1467      !   write(unitnb,*)    itau,opertau(itau)
1468      ! enddo
1469      ! write(unitnb,*) 
1470 
1471   ! -- Add correction for discontinuity.
1472 !      if(iflavor1==iflavor2) then
1473         !G(0+)-G(0-)=G(0+)+G(beta-)=A
1474         opertau(tauSamples+1) = -real(C) - opertau(1)
1475       !sui!write(6,*) "BackFourier",opertau(tauSamples+1),opertau(1),real(C)
1476 
1477         op%setT = .TRUE.
1478 !      endif
1479       op%oper(:,iflavor1,iflavor2)=opertau(:)
1480       FREE(opertau)
1481     END DO ! iflavor2
1482   END DO ! iflavor1
1483   ! -- End loop over flavors.
1484 
1485   FREE(Domega)
1486   FREE(A_omega)
1487   FREE(C_omega)
1488   close(236)
1489   close(237)
1490 
1491 END SUBROUTINE GreenHyboffdiag_backFourier

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_clear [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_clear

FUNCTION

  clear green function

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

388 SUBROUTINE GreenHyboffdiag_clear(op)
389 
390 !Arguments ------------------------------------
391 
392 !This section has been created automatically by the script Abilint (TD).
393 !Do not modify the following lines by hand.
394 #undef ABI_FUNC
395 #define ABI_FUNC 'GreenHyboffdiag_clear'
396 !End of the abilint section
397 
398   TYPE(GreenHyboffdiag)     , INTENT(INOUT) :: op
399   INTEGER :: iflavor,iflavorbis
400 
401   !CALL Vector_clear(op%oper_old)
402   !CALL VectorInt_clear(op%index_old)
403   do iflavor=1,op%nflavors
404     do iflavorbis=1,op%nflavors
405       CALL MapHyb_clear(op%map(iflavor,iflavorbis))
406     enddo
407   enddo
408   op%measurements = 0
409   IF ( ALLOCATED(op%oper) ) &
410   op%oper         = 0.d0
411   op%signvaluemeas = 0.d0
412   op%signvalueold = 1.d0
413   IF ( op%iTech .EQ. GREENHYB_OMEGA ) THEN
414     IF ( ALLOCATED(op%oper_w) ) &
415     op%oper_w       = CMPLX(0.d0,0.d0,8)
416     IF ( ALLOCATED(op%oper_w_old) ) &
417     op%oper_w_old   = CMPLX(0.d0,0.d0,8)
418   END IF
419   op%factor       = 0
420 END SUBROUTINE GreenHyboffdiag_clear

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_destroy [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_destroy

FUNCTION

  destroy green function

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

1928 SUBROUTINE GreenHyboffdiag_destroy(op)
1929 
1930 !Arguments ------------------------------------
1931 
1932 !This section has been created automatically by the script Abilint (TD).
1933 !Do not modify the following lines by hand.
1934 #undef ABI_FUNC
1935 #define ABI_FUNC 'GreenHyboffdiag_destroy'
1936 !End of the abilint section
1937 
1938   TYPE(GreenHyboffdiag), INTENT(INOUT) :: op
1939   INTEGER :: iflavor,iflavorbis
1940 
1941   op%set          = .FALSE.
1942   op%setT         = .FALSE.
1943   op%setW         = .FALSE.
1944   op%samples      = 0
1945   op%measurements = 0
1946   op%beta         = 0.d0
1947   op%inv_beta     = 0.d0
1948   op%inv_dt       = 0.d0
1949   op%delta_t      = 0.d0
1950   CALL VectorInt_destroy(op%index_old)
1951   CALL Vector_destroy(op%oper_old)
1952   do iflavor=1,op%nflavors
1953     do iflavorbis=1,op%nflavors
1954      !sui!write(6,*) "test",iflavor,iflavorbis
1955       CALL MapHyb_destroy(op%map(iflavor,iflavorbis))
1956     enddo
1957   enddo
1958   DT_FREEIF(op%map)
1959   FREEIF(op%oper)
1960   FREEIF(op%Mk)
1961   FREEIF(op%oper_w)
1962   FREEIF(op%oper_w_old)
1963   FREEIF(op%omega)
1964 END SUBROUTINE GreenHyboffdiag_destroy

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_forFourier [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_forFourier

FUNCTION

  perform forward fourier transform

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green
  Wmax=linear maximum frequency

OUTPUT

  Gomega=Results for omega frequencies
  omega=ask frequencies

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

1527 SUBROUTINE GreenHyboffdiag_forFourier(op, Gomega, omega, Wmax)
1528 !Arguments ------------------------------------
1529 
1530 !This section has been created automatically by the script Abilint (TD).
1531 !Do not modify the following lines by hand.
1532 #undef ABI_FUNC
1533 #define ABI_FUNC 'GreenHyboffdiag_forFourier'
1534 !End of the abilint section
1535 
1536 
1537 #ifdef HAVE_MPI1
1538 include 'mpif.h'
1539 #endif
1540   TYPE(GreenHyboffdiag)             , INTENT(INOUT) :: op
1541   COMPLEX(KIND=8), DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT) :: Gomega  ! INOUT for MPI
1542   COMPLEX(KIND=8), DIMENSION(:), OPTIONAL, INTENT(IN   ) :: omega  
1543   INTEGER                 , OPTIONAL, INTENT(IN   ) :: Wmax   
1544   INTEGER :: i
1545   INTEGER :: j
1546   INTEGER :: iflavor1
1547   INTEGER :: iflavor2
1548   INTEGER :: nflavors
1549   INTEGER :: L
1550   INTEGER :: Lspline
1551   INTEGER :: Nom
1552   INTEGER :: omegaBegin
1553   INTEGER :: omegaEnd
1554   INTEGER :: deltaw
1555   INTEGER :: residu
1556   INTEGER, ALLOCATABLE, DIMENSION(:) :: counts
1557   INTEGER, ALLOCATABLE, DIMENSION(:) :: displs
1558   DOUBLE PRECISION :: beta
1559   DOUBLE PRECISION :: tau
1560   DOUBLE PRECISION :: delta
1561   DOUBLE PRECISION :: deltabis
1562   DOUBLE PRECISION :: inv_delta
1563   DOUBLE PRECISION :: inv_delta2
1564   DOUBLE PRECISION :: omdeltabis
1565   DOUBLE PRECISION :: tmp
1566   DOUBLE PRECISION :: xpi
1567   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  diag
1568   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  diagL
1569   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  lastR
1570   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  lastC
1571   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  XM
1572   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  X2
1573   DOUBLE PRECISION :: iw
1574   COMPLEX(KIND=8) :: iwtau
1575   COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:) :: Gwtmp  
1576   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: omegatmp
1577   nflavors=op%nflavors
1578 
1579 !sui!write(6,*) " Fourier transformation begin"
1580 
1581   IF ( op%set .EQV. .FALSE. ) &
1582     CALL ERROR("GreenHyboffdiag_forFourier : Uninitialized GreenHyboffdiag structure")
1583   IF ( op%setT .EQV. .FALSE. ) &
1584     CALL ERROR("GreenHyboffdiag_forFourier : no G(tau)")
1585   !write(6,*) "op%setMk=", op%setMk
1586   IF ( op%setMk .NE. 2*nflavors*nflavors ) &
1587     CALL WARNALL("GreenHyboffdiag_forFourier : green does not have moments    ")
1588 
1589   L  = op%samples
1590 
1591   xpi=acos(-1.d0)                !!! XPI=PI
1592   beta = op%beta
1593   Nom  = op%Wmax
1594   IF ( PRESENT(Gomega) ) THEN
1595     Nom = SIZE(Gomega,1)    
1596     !IF ( op%rank .EQ. 0 ) &
1597       !!write(6,*) "size Gomega", Nom
1598   END IF
1599   IF ( PRESENT(omega) ) THEN
1600     IF ( PRESENT(Gomega) .AND. SIZE(omega) .NE. Nom ) THEN
1601       CALL ERROR("GreenHyboffdiag_forFourier : sizes mismatch              ")               
1602     !ELSE 
1603       !Nom = SIZE(omega)
1604     END IF 
1605   END IF
1606   IF ( .NOT. PRESENT(Gomega) .AND. .NOT. PRESENT(omega) ) THEN
1607     IF ( PRESENT(Wmax) ) THEN
1608       Nom=Wmax
1609     ELSE
1610       CALL ERROR("GreenHyboffdiag_forFourier : Missing argument Wmax")
1611     END IF
1612   END IF
1613 
1614   !!IF ( ALLOCATED(op%oper_w) ) THEN
1615   !!  IF ( SIZE(op%oper_w,1) .NE. Nom ) THEN
1616   !!    FREE(op%oper_w)
1617   !!    MALLOC(op%oper_w,(1:Nom,nflavors,nflavors))
1618   !!  END IF
1619   !!ELSE
1620   !!  MALLOC(op%oper_w,(1:Nom,nflavors,nflavors))
1621   !!END IF
1622 
1623   !!write(6,*) "PRESENT(GOMEGA)", PRESENT(GOMEGA)
1624   !!write(6,*) "PRESENT(OMEGA)", PRESENT(OMEGA)
1625   !call flush(6)
1626 
1627   delta=op%delta_t
1628   inv_delta = op%inv_dt
1629   inv_delta2 = inv_delta*inv_delta
1630  
1631   MALLOC(diagL,(L-1))
1632   MALLOC(lastR,(L-1))
1633   MALLOC(diag,(L))
1634   MALLOC(lastC,(L-1))
1635 
1636 !(cf Stoer) for the spline interpolation : 
1637 ! second derivatives XM solution of A*XM=B.
1638 !A=(2.4.2.11) of Stoer&Bulirsch + 2 limit conditions
1639 !The LU decomposition of A is known explicitly; 
1640 
1641   diag (1) = 4.d0 ! 1.d0 *4.d0 factor 4 added for conditionning
1642   diagL(1) = 0.25d0 !1.d0/4.d0
1643   lastR(1) = -0.5d0 ! -2.d0/4.d0
1644   lastC(1) = 4.d0 ! 1.d0*4.d0
1645   diag (2) = 4.d0
1646   diagL(2) = 0.25d0
1647   lastR(2) = -0.25d0
1648   lastC(2) = -1.d0
1649 
1650   DO i = 3, L-2
1651     tmp = 4.d0 - diagL(i-1)
1652     diagL(i) = 1.d0 / tmp
1653   END DO
1654   DO i = 3, L-2
1655     diag (i) = 1.d0 / diagL(i)
1656     lastR(i) = -(lastR(i-1)*diagL(i))
1657     lastC(i) = -(lastC(i-1)*diagL(i-1))
1658   END DO
1659   
1660   tmp = 1.d0/diag(L-2)
1661   diag (L-1) = 4.d0 - tmp
1662   lastR(L-1) = (1.d0 - lastR(L-2))/ diag(L-1)
1663   !diagL(L-1) = lastR(L-1)
1664   diagL(L-1) = 0.d0 ! for the Lq=B resolution
1665   !lastC(L-1) = 1.d0 - lastC(L-2)*diagL(L-1) ! equivalent to the next line
1666   lastC(L-1) = 1.d0 - (lastC(L-2)*lastR(L-1)) ! True value
1667   diag (L  ) = 2.d0! - DOT_PRODUCT( lastR , lastC )
1668   tmp = 0.d0
1669   DO i = 1, L-1
1670     tmp = tmp + lastR(i)*lastC(i)
1671   END DO
1672   diag (L  ) = diag (L  ) - tmp
1673   lastC(L-1) = lastC(L-1)-1.d0 ! 1 is removed for the u.XM=q resolution
1674 
1675   MALLOC(XM,(L))
1676   MALLOC(Gwtmp,(1:Nom))
1677 
1678   Lspline = L-1
1679   MALLOC(X2,(1:Lspline+1)) ! We impose L = Nom
1680 
1681   IF ( op%have_MPI .EQV. .TRUE. ) THEN  
1682     deltaw = Nom / op%size
1683     residu = Nom - op%size*deltaw
1684     IF ( op%rank .LT. op%size - residu ) THEN
1685       omegaBegin = 1 + op%rank*deltaw
1686       omegaEnd   = (op%rank + 1)*deltaw
1687     ELSE
1688   !    tauBegin = (op%size-residu)*deltaw + 1 + (op%rank-op%size+residu)*(deltaw+1)
1689       omegaBegin = 1 + op%rank*(deltaw + 1) -op%size + residu
1690       omegaEnd = omegaBegin + deltaw
1691     END IF
1692     MALLOC(counts,(1:op%size))
1693     MALLOC(displs,(1:op%size))
1694     counts = (/ (deltaw, i=1, op%size-residu), &
1695                 (deltaw+1, i=op%size-residu+1, op%size) /)
1696     displs(1)=0
1697     DO i = 2, op%size
1698       displs(i) = displs(i-1) + counts (i-1)
1699     END DO
1700   ELSE
1701     omegaBegin = 1
1702     omegaEnd   = Nom 
1703   END IF
1704 
1705 !  op%Mk(iflavor1,iflavor2,1) = 0.d0
1706 !  DO iflavor1 = 1, nflavors
1707 !    op%Mk(iflavor1,iflavor1,1) = -1.d0
1708 !  ENDDO
1709   op%Mk(:,:,3) = 0.d0
1710 
1711   MALLOC(omegatmp,(omegaBegin:omegaEnd))
1712   IF ( PRESENT(omega) ) THEN
1713     omegatmp(omegaBegin:omegaEnd) = (/ (AIMAG(omega(i)),i=omegaBegin,omegaEnd) /)
1714   ELSE
1715     omegatmp(omegaBegin:omegaEnd) = (/ ((((2.d0*DBLE(i)-1.d0)*xpi)/Beta), i=omegaBegin,omegaEnd) /)
1716   END IF
1717 
1718   DO iflavor1 = 1, nflavors
1719     DO iflavor2 = 1, nflavors
1720    ! write(6,*) "   Moments:",op%Mk(iflavor1,iflavor2,:),iflavor1,iflavor2
1721   
1722 ! construct the B vector from A.Xm=B
1723       XM(1) = 4.d0*op%Mk(iflavor1,iflavor2,3)
1724       XM(L) = (6.d0 * inv_delta) * ( op%Mk(iflavor1,iflavor2,2) - ( &
1725         (op%oper(2,iflavor1,iflavor2)-op%oper(1,iflavor1,iflavor2)) + &
1726         (op%oper(L,iflavor1,iflavor2)-op%oper(L-1,iflavor1,iflavor2)) ) * inv_delta )
1727 !    built generic second derivative of oper
1728 !sui!write(6,*)  "XM 1 L",XM(1),XM(L),op%Mk(iflavor1,iflavor2,2),op%Mk(iflavor1,iflavor2,3)
1729       DO i = 2, L-1
1730         XM(i) = (6.d0 * inv_delta2) * ( (op%oper(i+1,iflavor1,iflavor2) &
1731           - 2.d0 * op%oper(i,iflavor1,iflavor2)) &
1732           +        op%oper(i-1,iflavor1,iflavor2) )
1733     !sui!write(6,*) "XM",i,XM(i),op%oper(i,iflavor1,iflavor2)
1734       END DO
1735 
1736 ! Find second derivatives XM: Solve the system
1737 ! SOLVING Lq= XM 
1738 !  q = XM 
1739       do j=1,L-1
1740           XM(j+1)=XM(j+1)-(diagL(j)*XM(j))
1741           XM(L)  =XM(L)  -(lastR(j)*XM(j))
1742       end do
1743 
1744 
1745 ! SOLVING U.XM=q 
1746 !  XM = q
1747       do j=L-1,2,-1
1748        XM(j+1)  = XM(j+1) / diag(j+1)
1749        XM(j)= (XM(j)-(XM(L)*lastC(j)))-XM(j+1)
1750       end do
1751       XM(2)  = XM(2) / diag(2)
1752       XM(1) = (XM(1)-XM(L)*lastC(1)) / diag(1)
1753 
1754 
1755 
1756       !Construct L2 second derivative from known derivatives XM
1757       deltabis = beta / DBLE(Lspline)
1758       DO i = 1, Lspline
1759         tau = deltabis * DBLE(i-1)
1760         j = ((L-1)*(i-1))/Lspline + 1!INT(tau * inv_delta) + 1
1761         X2(i) = inv_delta * ( XM(j)*(DBLE(j)*delta - tau ) + XM(j+1)*(tau - DBLE(j-1)*delta) )
1762       END DO
1763       X2(Lspline+1) = XM(L)
1764 
1765 
1766        DO i = omegaBegin, omegaEnd
1767          iw = omegatmp(i)
1768          omdeltabis = iw*deltabis
1769          Gwtmp(i)=CMPLX(0.d0,0.d0,8)
1770          DO j=2, Lspline ! We impose  L+1 = Nom
1771            iwtau = CMPLX(0.d0,omdeltabis*DBLE(j-1),8)
1772            Gwtmp(i) = Gwtmp(i) + EXP(iwtau) * CMPLX((X2(j+1) + X2(j-1))-2.d0*X2(j),0.d0,8)
1773            !write(6,*) "ww",i,j,Gwtmp(i),X2(j),iwtau
1774          END DO
1775          Gwtmp(i) = Gwtmp(i)/CMPLX(((iw*iw)*(iw*iw)*deltabis),0.d0,8) &
1776          + CMPLX( ( ((X2(2)-X2(1))+(X2(Lspline+1)-X2(Lspline)))/((iw*iw)*deltabis) -op%Mk(iflavor1,iflavor2,2) ) &
1777          /(iw*iw) , (op%Mk(iflavor1,iflavor2,1)-op%Mk(iflavor1,iflavor2,3)/(iw*iw))/iw , 8) 
1778                    !+ CMPLX( (X2(2)-X2(1))+(X2(Lspline+1)-X2(Lspline)), 0.d0, 8 ) ) &
1779                    !   / (((iw*iw)*(iw*iw))*CMPLX(deltabis,0.d0,8)) &
1780                    !- CMPLX(op%Mk(1),0.d0,8)/iw  &
1781                    !+ CMPLX(op%Mk(2),0.d0,8)/(iw*iw) &
1782                    !- CMPLX(op%Mk(3),0.d0,8)/((iw*iw)*iw) 
1783          !IF ( op%rank .EQ. 0 )  write(12819,*) iw,gwtmp(i)
1784        END DO
1785        !call flush(12819)
1786        IF ( op%have_MPI .EQV. .TRUE. ) THEN
1787 #ifdef HAVE_MPI
1788         CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DATATYPE_NULL, &
1789                           Gwtmp  , counts, displs, &
1790                           MPI_DOUBLE_COMPLEX, op%MY_COMM, residu)
1791 #endif
1792       END IF
1793       IF ( PRESENT(Gomega) ) THEN
1794         Gomega(:,iflavor1,iflavor2) = Gwtmp(:)
1795       END IF
1796       op%setW = .TRUE.
1797     ENDDO ! iflavor1
1798   ENDDO ! iflavor2
1799   !!op%oper_w=Gomega
1800   do iflavor1=1,nflavors
1801     !sui!write(6,*)  iflavor1
1802       do i=1,Nom
1803       !write(6,*) "w",i,op%oper_w(i,iflavor1,iflavor1)
1804       enddo
1805   enddo
1806   FREE(Gwtmp)
1807   FREE(diagL)
1808   FREE(lastR)
1809   FREE(diag)
1810   FREE(lastC)
1811   FREE(XM)
1812   FREE(omegatmp)
1813   FREE(X2)
1814   FREE(counts)
1815   FREE(displs)
1816 
1817 END SUBROUTINE GreenHyboffdiag_forFourier

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_getHybrid [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_getHybrid

FUNCTION

  reduce green function

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

880 SUBROUTINE GreenHyboffdiag_getHybrid(op)
881 
882 !Arguments ------------------------------------
883 
884 !This section has been created automatically by the script Abilint (TD).
885 !Do not modify the following lines by hand.
886 #undef ABI_FUNC
887 #define ABI_FUNC 'GreenHyboffdiag_getHybrid'
888 !End of the abilint section
889 
890   TYPE(GreenHyboffdiag), INTENT(INOUT) :: op
891 
892   IF ( op%set .EQV. .FALSE. ) &
893     CALL ERROR("GreenHyboffdiag_getHybrid : green operator not set          ")
894 
895   SELECT CASE(op%iTech)
896   CASE (GREENHYB_TAU)
897     op%oper = -(op%oper * op%inv_beta) / (DBLE(op%measurements) * op%delta_t)
898   !sui!write(6,*) "measurements",op%measurements,op%delta_t,op%inv_beta
899   !sui!write(6,*) "signevaluemeas meas",op%signvaluemeas,op%measurements
900     op%signvaluemeas = op%signvaluemeas / DBLE(op%measurements)
901    ! print*, "op%oper",op%oper(1,1,1)
902   !sui!write(6,*) "signevaluemeas/meas",op%signvaluemeas
903    ! print*, "signevaluemeas/meas",op%signvaluemeas
904     op%setT = .TRUE.
905   CASE (GREENHYB_OMEGA)
906     op%oper_w = -(op%oper_w * op%inv_beta) / (DBLE(op%measurements) * op%delta_t)
907     op%setW = .TRUE.
908     CALL GreenHyboffdiag_backFourier(op)
909   END SELECT
910 
911 END SUBROUTINE GreenHyboffdiag_getHybrid

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_init [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_init

FUNCTION

  Initialize and allocate

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green
  samples=imaginary time slices
  beta=inverse temperature
  iTech=SHOULD NOT BE USED => BUGGY
  MY_COMM=mpi_communicator

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

201 SUBROUTINE GreenHyboffdiag_init(op, samples, beta,nflavors,iTech,MY_COMM)
202 
203 
204 !This section has been created automatically by the script Abilint (TD).
205 !Do not modify the following lines by hand.
206 #undef ABI_FUNC
207 #define ABI_FUNC 'GreenHyboffdiag_init'
208 !End of the abilint section
209 
210 
211 #ifdef HAVE_MPI1
212 include 'mpif.h'
213 #endif
214 !Arguments ------------------------------------
215   TYPE(GreenHyboffdiag)     , INTENT(INOUT) :: op
216   INTEGER         , INTENT(IN   ) :: samples
217   DOUBLE PRECISION, INTENT(IN   ) :: beta
218   INTEGER , INTENT(IN   ) :: nflavors
219   !INTEGER          , INTENT(IN   ) :: Wmax
220   INTEGER, OPTIONAL, INTENT(IN   ) :: iTech
221   INTEGER, OPTIONAL, INTENT(IN   ) :: MY_COMM
222 !Local variables ------------------------------
223   INTEGER                         :: iflavor,iflavorbis,sp1
224   DOUBLE PRECISION                :: dt
225 #ifdef HAVE_MPI
226   INTEGER          :: ierr
227 #endif
228 
229   IF ( PRESENT(MY_COMM)) THEN
230 #ifdef HAVE_MPI
231     op%have_MPI = .TRUE.
232     op%MY_COMM = MY_COMM
233     CALL MPI_Comm_rank(op%MY_COMM, op%rank, ierr)
234     CALL MPI_Comm_size(op%MY_COMM, op%size, ierr)
235 #else
236     CALL WARN("GreenHyboffdiag_init : MPI is not used                                    ")
237     op%have_MPI = .FALSE.
238     op%MY_COMM = -1
239     op%rank = 0
240     op%size = 1
241 #endif
242   ELSE
243     op%have_MPI = .FALSE.
244     op%MY_COMM = -1
245     op%rank = 0
246     op%size = 1
247   END IF
248 
249   sp1             = samples + 1
250   op%samples      = sp1
251   op%measurements = 0
252   op%nflavors     = nflavors
253   op%beta         = beta
254   op%inv_beta     = 1.d0 / beta
255   op%inv_dt       = DBLE(samples) * op%inv_beta
256   dt              = 1.d0 / op%inv_dt
257   op%delta_t      = dt
258   !op%Wmax         = Wmax
259   op%Wmax         = -1
260   FREEIF(op%Mk)
261   MALLOC(op%Mk,(nflavors,nflavors,3))
262   FREEIF(op%oper)
263   MALLOC(op%oper,(sp1,nflavors,nflavors))
264   ! If we want to measure in frequences
265   ! let assume we first have "samples" frequences
266   IF ( PRESENT(iTech) ) THEN
267     op%iTech = iTech
268     SELECT CASE (op%iTech)
269     CASE (GREENHYB_TAU)  ! omega
270       op%iTech = GREENHYB_TAU
271     CASE (GREENHYB_OMEGA)  ! omega
272       op%Wmax = samples
273       FREEIF(op%oper_w)
274       MALLOC(op%oper_w,(1:op%Wmax,nflavors,nflavors))
275       FREEIF(op%oper_w_old)
276       MALLOC(op%oper_w_old,(1:op%Wmax))
277       op%oper_w     = CMPLX(0.d0,0.d0,8)
278       op%oper_w_old = CMPLX(0.d0,0.d0,8)
279       FREEIF(op%omega)
280       MALLOC(op%omega,(1:op%Wmax))
281       op%omega = (/ ((2.d0 * DBLE(sp1) - 1.d0)*ACOS(-1.d0)*op%inv_beta, sp1=1, op%Wmax) /)
282     END SELECT
283   ELSE
284     op%iTech = GREENHYB_TAU
285   END IF
286   ! end if
287   CALL Vector_init(op%oper_old,10000)
288   CALL VectorInt_init(op%index_old,10000)
289   DT_FREEIF(op%map)
290   MALLOC(op%map,(nflavors,nflavors))
291   do iflavor=1,nflavors
292     do iflavorbis=1,nflavors
293       CALL MapHyb_init(op%map(iflavor,iflavorbis),10000)
294     enddo
295   enddo
296 
297   op%oper       = 0.d0
298   op%signvaluemeas = 0.d0
299   op%signvalueold = 0.d0
300   op%set        = .TRUE.
301   op%factor     = 1
302   op%setMk      = 0
303   op%Mk         = 0.d0
304 END SUBROUTINE GreenHyboffdiag_init

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_measHybrid [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_measHybrid

FUNCTION

  Measure Green's function

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green
  Mmatrix=M matrix for the current flavor
  ListCdagC_1=list of all creator and annhilator operators
  updated=should we accumulate or not

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

519 SUBROUTINE GreenHyboffdiag_measHybrid(op, Mmatrix, ListCdagC_1, updated,signvalue,activeflavor)
520 
521 !Arguments ------------------------------------
522 
523 !This section has been created automatically by the script Abilint (TD).
524 !Do not modify the following lines by hand.
525 #undef ABI_FUNC
526 #define ABI_FUNC 'GreenHyboffdiag_measHybrid'
527 !End of the abilint section
528 
529   TYPE(GreenHyboffdiag)    , INTENT(INOUT) :: op
530   TYPE(MatrixHyb)   , INTENT(IN   ) :: Mmatrix
531   TYPE(ListCdagC)   , INTENT(IN   ) :: ListCdagC_1(op%nflavors)
532   DOUBLE PRECISION  , INTENT(IN   ) :: signvalue
533   LOGICAL        , INTENT(IN   ) :: updated
534   INTEGER, OPTIONAL  , INTENT(IN  ) :: activeflavor
535 !Local variables ------------------------------
536   INTEGER                        :: iC
537   INTEGER                        :: iCdag
538   INTEGER                        :: tail
539   INTEGER                        :: tailbis
540   !INTEGER                        :: index
541   INTEGER                        :: idx_old
542   INTEGER                        :: old_size
543 !  INTEGER                        :: omegaSamples
544 !  INTEGER                        :: iomega
545   INTEGER                        :: iflavor
546   INTEGER                        :: iflavorbis
547   INTEGER                        :: iC_m,iC_m_add
548   INTEGER                        :: iCdag_m,iCdag_m_add
549   INTEGER                        :: stail,ii
550 !  DOUBLE PRECISION               :: pi_invBeta
551   DOUBLE PRECISION               :: mbeta_two
552   DOUBLE PRECISION               :: beta 
553   DOUBLE PRECISION               :: beta_tc
554   DOUBLE PRECISION               :: tcbeta_tc
555   DOUBLE PRECISION               :: inv_dt 
556   DOUBLE PRECISION               :: tC,tc_phys
557   DOUBLE PRECISION               :: tCdag
558   DOUBLE PRECISION               :: time
559   DOUBLE PRECISION               :: signe,signe2
560   DOUBLE PRECISION               :: argument
561   INTEGER                        :: iflavorbegin,iflavorend,prtopt
562   !DOUBLE PRECISION               :: taupi_invbeta
563   !COMPLEX(KIND=8)                   :: cargument
564   !COMPLEX(2*8)                   :: base_exp
565   !COMPLEX(2*8)                   :: increm_exp
566   !write(6,*) "measHybrid"
567   prtopt=0
568   IF ( op%set .EQV. .FALSE. ) &
569     CALL ERROR("GreenHyboffdiag_measHybrid : green operator not set         ")
570   stail=0
571   do iflavor=1,op%nflavors
572     stail=stail + ListCdagC_1(iflavor)%tail
573   enddo
574   iflavorbegin = 1
575   iflavorend   = op%nflavors
576 
577   if(present(activeflavor)) then
578     if(activeflavor.ne.0) then
579 !sui!write(6,*) "measHybrid activeflavor",activeflavor
580       iflavorbegin = activeflavor
581       iflavorend   = activeflavor
582     endif
583   endif
584 
585   IF ( stail .NE. Mmatrix%tail ) &
586       CALL ERROR("GreenHyboffdiag_measHybrid : ListCdagC & M unconsistent     ")
587 
588   IF ( updated .EQV. .TRUE. ) THEN ! NEW change in the configuration
589     ! FIXME SHOULD be much more faster
590     
591 
592     SELECT CASE(op%iTech)
593     CASE (GREENHYB_TAU)
594       argument = DBLE(op%factor)
595 !     At the beginning old_size=0, then it increases
596 !     until 
597 !     for all values of iC, increment green%oper with the value of the
598 !     Green's function in listDBLE(iC) obtained from previous iteration
599 !     (below)
600       ! ===============================================================
601       ! An update has been done. So the Green's function will change
602       ! It is thus the good moment to store the previous Green's
603       ! function with argument, the number of times this Green's
604       ! function has been constant
605       ! ===============================================================
606       DO iflavor=1, op%nflavors
607         DO iflavorbis=1, op%nflavors
608           old_size = op%map(iflavor,iflavorbis)%tail
609           !write(6,*) "size listDBLE",size(op%map(iflavor,iflavorbis)%listDBLE)
610         !sui!write(6,*) " measHybrid",old_size,iflavor,iflavorbis
611           DO iC = 1, old_size
612             if(op%map(iflavor,iflavorbis)%listINT(iC)==0) then
613               !write(6,*) "listINT(iC)=",iC,op%map(iflavor,iflavorbis)%listINT(iC)
614             endif
615             !write(6,*) " measHybrid  iflavor,iflavorbis,iC listINT ",iflavor,iflavorbis,iC,op%map(iflavor,iflavorbis)%listINT(iC)
616                !write(6,*) " measHybrid  listDBLE ",iflavor,iflavorbis,iC,op%map(iflavor,iflavorbis)%listDBLE(iC),argument
617             op%oper(op%map(iflavor,iflavorbis)%listINT(iC),iflavor,iflavorbis) =                &
618                            op%oper(op%map(iflavor,iflavorbis)%listINT(iC),iflavor,iflavorbis) &
619                          + op%map(iflavor,iflavorbis)%listDBLE(iC) * op%signvalueold *  argument
620            !if(op%map(iflavor,iflavorbis)%listINT(iC)==1.and.iflavor==iflavorbis) then
621           !  if(iflavor==iflavorbis) then
622           !   !sui!write(6,*) "G(0)", op%map(iflavor,iflavorbis)%listINT(iC),op%map(iflavor,iflavorbis)%listDBLE(iC) * op%signvalueold,op%oper(op%map(iflavor,iflavorbis)%listINT(iC),iflavor,iflavorbis),iflavor
623           !  endif
624           !  if(iflavor==1.and.iflavorbis==6.and.op%map(iflavor,iflavorbis)%listINT(iC)==1) then
625           !          !prt!if(prtopt==1) write(6,*) "G16(0)", op%map(iflavor,iflavorbis)%listDBLE(iC) * op%signvalueold*argument,op%oper(op%map(iflavor,iflavorbis)%listINT(iC),iflavor,iflavorbis),op%signvalueold
626           !  endif
627           !  if(iflavor==6.and.iflavorbis==1.and.op%map(iflavor,iflavorbis)%listINT(iC)==1) then
628           !          !prt!if(prtopt==1) write(6,*) "G61(0)", op%map(iflavor,iflavorbis)%listDBLE(iC) * op%signvalueold*argument,op%oper(op%map(iflavor,iflavorbis)%listINT(iC),iflavor,iflavorbis),op%signvalueold
629           !  endif
630           END DO
631       ! tail**2 is the number of possible t-t'
632       ! MapHyb_setSize with resize map tail*tail will thus be the new
633       ! op%map%tail
634       ! update size of map and map%tail
635           CALL MapHyb_setSize(op%map(iflavor,iflavorbis),&
636 &          ListCdagC_1(iflavor)%tail*ListCdagC_1(iflavorbis)%tail)
637         END DO
638       END DO
639       op%signvaluemeas = op%signvaluemeas + op%signvalueold * argument
640       op%measurements = op%measurements + op%factor
641     !sui!write(6,*) "   measurements", op%measurements
642          !sui! write(6,*) "                  signvaluemeas",op%signvaluemeas,op%signvalueold*argument
643          !sui! write(6,*) "                  signvaluemeas/measurements",op%signvaluemeas/op%measurements
644   
645       ! This is new measurement, thus op%factor should be put to one
646       op%factor = 1
647 
648 
649       ! initialized index idx_old for the doubles loops over flavors and segments.
650 
651       ! setup usefull quantities
652       beta   =  op%beta
653       mbeta_two = -(beta*0.5d0)
654       inv_dt =  op%inv_dt
655 
656       ! WARNING time is not the time but just a temporary variable.
657       ! Index Time has been calculated previously and is in mat_tau
658 
659       ! initialized index for each annihilation time of a segment for a given flavor
660 
661       ! initialized index for each creation time of a segment for another  flavor
662 
663       iC_m=0
664       iC_m_add=0
665       DO iflavor=1,op%nflavors
666         tail=ListCdagC_1(iflavor)%tail
667         !write(6,*) " measHybrid  iflavor",iflavor,tail
668 
669         iCdag_m=0
670         iCdag_m_add=0
671         DO iflavorbis=1,op%nflavors
672           tailbis=ListCdagC_1(iflavorbis)%tail
673           !write(6,*) " measHybrid  iflavorbis",iflavorbis,tailbis
674           idx_old = 0
675 
676           DO iC  = 1, tail
677             ! tC is the annihilation (C_) time for segment iC and flavor iflavor
678             !-------------------------------------------------------------------
679             tC   = ListCdagC_1(iflavor)%list(iC,C_)
680 
681             !iC_m=iC_m+1  ! for Mmatrix%mat
682             ! For each flavor iflavor, iC start at \sum_{iflavor1<iflavor} tail(iflavor1)
683             ! It thus explains the presence of iC_m_add (same below for iCdag_m_add)
684             ! ---------------------------------------------------------------------------------
685             iC_m=iC_m_add+iC
686             beta_tc = beta - tC
687             tcbeta_tc = tC * beta_tc
688 
689              !write(6,*) " measHybrid  iC_m",iC_m
690              !write(6,*) " measHybrid  tailbis",tailbis
691             DO iCdag = 1, tailbis
692               !iCdag_m=iCdag_m+1
693               iCdag_m=iCdag_m_add+iCdag
694              !write(6,*) " measHybrid  iCdag_m",iCdag_m
695 
696               ! tCdag is the creation time for segment iCdag and flavor iflavorbis
697               tCdag  = ListCdagC_1(iflavorbis)%list(iCdag,Cdag_)
698 
699 !  ---        time is equivalent to  time=(tc-tcdag)*(beta-tc) and is only
700 !  ---        useful for signe
701               time = tcbeta_tc - tCdag*beta_tc
702   
703               !signe = SIGN(1.d0,time)
704               !time = time + (signe-1.d0)*mbeta_two
705               !signe = signe * SIGN(1.d0,beta-tC)
706               !signe = SIGN(1.d0,time) * SIGN(1.d0,beta-tC)
707               tc_phys=tc
708               if(tc>beta) tc_phys=tc-beta
709               signe2=SIGN(1.d0,tc_phys-tcdag)
710 
711               if(iflavor==iflavorbis) signe = SIGN(1.d0,time)
712               if(iflavor/=iflavorbis) signe = signe2
713              ! signe = SIGN(1.d0,tc-tcdag)
714               ! --- tc>tcdag and beta>tc signe=1  ! segment in the middle or antisegment at the edge
715               !                                   ! tc-tcdag  > 0
716               ! --- tc<tcdag and beta<tc signe=1  ! never
717               ! --- tc>tcdag and beta<tc signe=-1 ! segment  at the edges
718               !                                   ! tc'-tcdag < 0 (with tc'=tc-beta)  -> signe < 0
719               ! --- tc<tcdag and beta>tc signe=-1 ! antisegment in the middle 
720               !                                   ! tc-tcdag  < 0 (with tc'=tc-beta)  -> signe < 0
721               ! 22/09/14:
722               ! ListCdagC_1 is the list of segment, so we are dealing
723               ! only with segment here. However all combination of Cdag
724               ! and C are taken, this it is possible that tc<tcdag
725 
726               ! 21/10/14: Wagt is important are the true times (between
727               ! 0 and beta). If tauC>tauCdag signe=+1
728               !              If tauC<tauCdag signe=-1
729               ! if(tc<tcdag.and.(iflavor==iflavorbis)) then
730               !   write(6,*)  ListCdagC_1(iflavorbis)%tail
731               !   do ii=1, ListCdagC_1(iflavorbis)%tail
732               !     write(6,*)  ListCdagC_1(iflavorbis)%list(ii,1), ListCdagC_1(iflavorbis)%list(ii,2)
733               !   enddo
734               !   write(6,*) "tc<tcdag", tc,tcdag,beta,iflavor,iflavorbis
735               !   stop
736               ! endif
737 
738               if(tc-tcdag>beta) then
739              !   write(6,*) " tc-tcdag > beta ", tcdag-tc,beta
740               endif
741               !if(tc>beta) then
742               !  write(6,*) " TC>BETA"
743               !  write(6,*) " iflavor,iflavorbis",iflavor,iflavorbis
744               !  write(6,*) " ic,icdag          ",ic,icdag
745               !  write(6,*) " signe             ",signe
746               !  write(6,*) " tc,tcdag          ",tc,tcdag
747               !  write(6,*) " Mmatrix%mat       ",Mmatrix%mat(iCdag_m,iC_m)
748               !  write(6,*) " Mmatrix%mat_tau   ",Mmatrix%mat_tau(iCdag_m,iC_m)
749               !endif
750               ! Si iflavor/=iflavorbis, tc-tcdag can be negative..so in
751               ! this case, on should add beta to tc-tcdag with the minus
752               ! sign. NOT DONE HERE??
753   
754 !             ----- Compute the Green's function as the value of the matrix M for times iCdag and iC.
755               argument = signe*Mmatrix%mat(iCdag_m,iC_m)
756   
757               !index = INT( ( time * inv_dt ) + 1.5d0 )
758               !IF (index .NE. Mmatrix%mat_tau(iCdag,iC)) THEN
759               !  WRITE(*,*) index, Mmatrix%mat_tau(iCdag,iC)
760               !!  CALL ERROR("Plantage")
761               !END IF
762   
763               idx_old = idx_old + 1
764 
765               ! --- define the  value of listDBLE as a function of idx_old
766               op%map(iflavor,iflavorbis)%listDBLE(idx_old) = argument
767               !write(6,*) " measHybrid  listDBLE2 ",iflavor,iflavorbis,idx_old,argument
768               !op%map%listINT(idx_old)  = index
769 
770               ! --- define the new corresponding value of listINT(idx_old) from mat_tau (integers)
771               ! --- idx_old has no meaning but listINT(idx_old) has.
772               op%map(iflavor,iflavorbis)%listINT(idx_old)  = Mmatrix%mat_tau(iCdag_m,iC_m)
773              !write(6,*) " measHybrid  idx_old listINT ",idx_old,op%map(iflavor,iflavorbis)%listINT(idx_old)
774              !write(6,*) " measHybrid  iCdag_m, iC_m, mat_tau",iCdag_m,iC_m,Mmatrix%mat_tau(iCdag_m,iC_m)
775             !  if(iflavor==1.and.iflavorbis==2.and.op%map(iflavor,iflavorbis)%listINT(idx_old)==1) then
776             !    !prt!if(prtopt==1) write(6,*) "---------------------------"
777             !    !prt!if(prtopt==1) write(6,*) "GG12(0)", op%map(iflavor,iflavorbis)%listINT(idx_old),op%map(iflavor,iflavorbis)%listDBLE(idx_old),tcdag,tc,signe,signe2
778             !    !prt!if(prtopt==1) write(6,*) "       ", tc-tcdag,tc_phys-tcdag
779             !    do ii=1, tail
780             !      !prt!if(prtopt==1)  write(6,*) ii, ListCdagC_1(iflavor)%list(ii,1), ListCdagC_1(iflavor)%list(ii,2)
781             !    enddo
782             !    do ii=1, tailbis
783             !      !prt!if(prtopt==1)  write(6,*) ii, ListCdagC_1(iflavorbis)%list(ii,1), ListCdagC_1(iflavorbis)%list(ii,2)
784             !    enddo
785             !    !prt!if(prtopt==1) write(6,*) "---------------------------"
786             !  endif
787             !  if(iflavor==1.and.iflavorbis==2.and.op%map(iflavor,iflavorbis)%listINT(idx_old)==999) then
788             !    !prt!if(prtopt==1) write(66,*) "---------------------------"
789             !    !prt!if(prtopt==1) write(66,*) "GG12(0)", op%map(iflavor,iflavorbis)%listINT(idx_old),op%map(iflavor,iflavorbis)%listDBLE(idx_old),tcdag,tc,signe,signe2
790             !    !prt!if(prtopt==1) write(66,*) "       ", tc-tcdag,tc_phys-tcdag
791             !    do ii=1, tail
792             !      !prt!if(prtopt==1)  write(66,*) ii, ListCdagC_1(iflavor)%list(ii,1), ListCdagC_1(iflavor)%list(ii,2)
793             !    enddo
794             !    do ii=1, tailbis
795             !      !prt!if(prtopt==1)  write(66,*) ii, ListCdagC_1(iflavorbis)%list(ii,1), ListCdagC_1(iflavorbis)%list(ii,2)
796             !    enddo
797             !    !prt!if(prtopt==1) write(66,*) "---------------------------"
798             !  endif
799             !  !if(iflavor==2.and.iflavorbis==1.and.op%map(iflavor,iflavorbis)%listINT(idx_old)==1) then
800               !   !prt!if(prtopt==1) write(6,*) "GG21(0)", op%map(iflavor,iflavorbis)%listINT(idx_old),op%map(iflavor,iflavorbis)%listDBLE(idx_old),tcdag,tc,signe
801               !endif
802 
803 
804             END DO
805           END DO
806         !  do ii=1,tail*tailbis
807         !   !write(6,*) " measHybrid  ii,op%map(iflavor,iflavorbis)%listINT(ii)", ii,op%map(iflavor,iflavorbis)%listINT(ii)
808         !  enddo
809           iCdag_m_add=iCdag_m_add+tailbis
810         END DO ! iflavorbis
811        iC_m_add=iC_m_add+tail
812       END DO ! iflavor
813       op%signvalueold = signvalue
814     CASE (GREENHYB_OMEGA)
815     !  argument = DBLE(op%factor)
816     !  DO iomega = 1, omegaSamples
817     !    op%oper_w(iomega) = op%oper_w(iomega) + op%oper_w_old(iomega) * argument
818     !  END DO
819     !  op%measurements = op%measurements + op%factor
820 
821     !  op%factor = 1
822     !  beta   =  op%beta
823     !  mbeta_two = -(beta*0.5d0)
824     !  pi_invBeta = ACOS(-1.d0)/beta
825     !  omegaSamples = op%samples-1
826     !  DO iC  = 1, tail
827     !    tC   = ListCdagC_1%list(iC,C_)
828     !    DO iCdag = 1, tail
829     !      tCdag  = ListCdagC_1%list(iCdag,Cdag_)
830     !      time = tC - tCdag
831 
832     !      signe = SIGN(1.d0,time)
833     !      time = time + (signe-1.d0)*mbeta_two
834     !      signe = signe * SIGN(1.d0,beta-tC)
835     !      argument = signe*Mmatrix%mat(iCdag,iC)
836 
837     !      DO iomega = 1, omegaSamples
838     !        !op%oper_w_old(iomega) = Mmatrix%mat_tau(iCdag,iC)*CMPLX(0.d0,argument)
839     !        op%oper_w_old(iomega) = EXP(CMPLX(0.d0,op%omega(iomega)*time))*CMPLX(0.d0,argument)
840     !      END DO
841     !    END DO
842     !  END DO
843     END SELECT
844   ELSE
845     op%factor = op%factor + 1
846   END IF
847 END SUBROUTINE GreenHyboffdiag_measHybrid

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_print [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_print

FUNCTION

  print Green function

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green
  ostream=file stream

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

1851 SUBROUTINE GreenHyboffdiag_print(op, ostream)
1852 
1853 !Arguments ------------------------------------
1854 
1855 !This section has been created automatically by the script Abilint (TD).
1856 !Do not modify the following lines by hand.
1857 #undef ABI_FUNC
1858 #define ABI_FUNC 'GreenHyboffdiag_print'
1859 !End of the abilint section
1860 
1861   TYPE(GreenHyboffdiag), INTENT(IN) :: op
1862   INTEGER, OPTIONAL , INTENT(IN) :: ostream
1863 !Local variables ------------------------------
1864   INTEGER                        :: ostream_val
1865   INTEGER                        :: isample
1866   INTEGER                        :: samples
1867   INTEGER                        :: iflavor1
1868   INTEGER                        :: iflavor2
1869   
1870 
1871   IF ( op%set .EQV. .FALSE. ) &
1872     CALL ERROR("GreenHyboffdiag_print : green op%operator not set              ")
1873 
1874   IF ( PRESENT(ostream) ) THEN
1875     ostream_val = ostream
1876   ELSE
1877     ostream_val = 66
1878     OPEN(UNIT=ostream_val,FILE="Green.dat")
1879   END IF
1880 
1881   samples =  op%samples 
1882 
1883   DO iflavor1=1,op%nflavors
1884     DO iflavor2=1,op%nflavors
1885     WRITE(ostream_val,'(a,i3,a,i3,a)')  "## (iflavor1,iflavor2)= (", iflavor1,",",iflavor2,")"
1886       DO isample = 1, samples
1887       WRITE(ostream_val,*) DBLE(isample-1)*op%delta_t, op%oper(isample,iflavor1,iflavor2)
1888       END DO
1889       WRITE(ostream_val,*) 
1890     END DO
1891   END DO
1892 
1893   IF ( .NOT. PRESENT(ostream) ) &
1894     CLOSE(ostream_val)
1895 END SUBROUTINE GreenHyboffdiag_print

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_reset [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_reset

FUNCTION

  reset green function

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

337 SUBROUTINE GreenHyboffdiag_reset(op)
338 
339 !Arguments ------------------------------------
340 
341 !This section has been created automatically by the script Abilint (TD).
342 !Do not modify the following lines by hand.
343 #undef ABI_FUNC
344 #define ABI_FUNC 'GreenHyboffdiag_reset'
345 !End of the abilint section
346 
347   TYPE(GreenHyboffdiag)     , INTENT(INOUT) :: op
348 
349   CALL GreenHyboffdiag_clear(op)
350   op%setMk        = 0
351   op%Mk           = 0.d0
352   op%setT         = .FALSE.
353   op%setW         = .FALSE.
354 END SUBROUTINE GreenHyboffdiag_reset

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setMoments [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_setMoments

FUNCTION

  Compute full moments

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Greenb
  u1_iflavor1=\sum_{iflavor2} U_{iflavor2,iflavor1} N_iflavor2
    (useful for first moment)  
  u2=\sum_{iflavor1,iflavor2,iflavor3} U_{iflavor1,iflavor2} N_iflavor2

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

 CHI
  Will be filled automatically by the parent script

SOURCE

1096 SUBROUTINE GreenHyboffdiag_setMoments(op,iflavor1,iflavor1b,u1,u2,u3)
1097 
1098 !Arguments ------------------------------------
1099 
1100 !This section has been created automatically by the script Abilint (TD).
1101 !Do not modify the following lines by hand.
1102 #undef ABI_FUNC
1103 #define ABI_FUNC 'GreenHyboffdiag_setMoments'
1104 !End of the abilint section
1105 
1106   TYPE(GreenHyboffdiag)  , INTENT(INOUT) :: op
1107   DOUBLE PRECISION, INTENT(IN   ) :: u1
1108   DOUBLE PRECISION, INTENT(IN   ) :: u2
1109   DOUBLE PRECISION, INTENT(IN   ) :: u3
1110   INTEGER         , INTENT(IN   ) :: iflavor1
1111   INTEGER         , INTENT(IN   ) :: iflavor1b
1112   
1113   if(iflavor1==iflavor1b) then
1114     op%Mk(iflavor1,iflavor1b,1) = -1.d0
1115 !   c_a(3)=-d1-mu*mu-2(-mu)(\sum_{b.ne.a} Uab nb) 
1116     op%Mk(iflavor1,iflavor1b,3) = op%Mk(iflavor1,iflavor1b,3) - 2.d0*(op%Mk(iflavor1,iflavor1b,2)*u1)
1117 
1118 !   c_a(2)=-mu+\sum_{b.ne.a} Uab n_b
1119     op%Mk(iflavor1,iflavor1b,2) = op%Mk(iflavor1,iflavor1b,2) + u1
1120  !sui!write(6,*) "setmiments",iflavor1,iflavor1b,u1
1121 
1122 !   c_a(3)=c_a(3) + \sum Uab^2 nb + \sum Uba Uca <nbnc> 
1123 !   ie c_a(3)=-d1+mu*mu-2mu*\sumb Uab nb + \sum Uab^2 nb + \sum Uba Uca <nbnc> 
1124     op%Mk(iflavor1,iflavor1b,3) = op%Mk(iflavor1,iflavor1b,3) - u2
1125   else
1126     op%Mk(iflavor1,iflavor1b,1) = 0.d0
1127     op%Mk(iflavor1,iflavor1b,2) = u3
1128     op%Mk(iflavor1,iflavor1b,3) = 0.d0
1129   endif
1130 !write(6,*) "mom",iflavor1,iflavor1b, op%Mk(iflavor1,iflavor1b,:)
1131 
1132   op%setMk = op%setMk + 1
1133 
1134 END SUBROUTINE GreenHyboffdiag_setMoments

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setMuD1 [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_setMuD1

FUNCTION

  Set first moments for G

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green
  mu=energy level (irrespectige with fermi level)
  d1=first moment of hybridization function ("K")

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

1024 SUBROUTINE GreenHyboffdiag_setMuD1(op,iflavor,iflavor2,mu,d1)
1025 
1026 !Arguments ------------------------------------
1027 !Arguments ------------------------------------
1028 !scalars
1029 
1030 !This section has been created automatically by the script Abilint (TD).
1031 !Do not modify the following lines by hand.
1032 #undef ABI_FUNC
1033 #define ABI_FUNC 'GreenHyboffdiag_setMuD1'
1034 !End of the abilint section
1035 
1036   DOUBLE PRECISION, INTENT(IN   ) :: mu
1037   DOUBLE PRECISION, INTENT(IN   ) :: d1
1038   INTEGER         , INTENT(IN   ) :: iflavor
1039   INTEGER         , INTENT(IN   ) :: iflavor2
1040 !type
1041   TYPE(GreenHyboffdiag)  , INTENT(INOUT) :: op
1042 !Local variables ------------------------------------
1043   DOUBLE PRECISION                :: mu2
1044 !*********************************************************************
1045   
1046   mu2=0
1047   if(iflavor==iflavor2) mu2=mu
1048 
1049   if(iflavor==iflavor2) then
1050     !op%Mk(iflavor,iflavor2,3) = -d1-(mu*mu)
1051     op%Mk(iflavor,iflavor2,3) = -(mu*mu)
1052     op%Mk(iflavor,iflavor2,2) = -mu
1053  !sui!write(6,*) "setmud1",iflavor,iflavor2, op%Mk(iflavor,iflavor2,2), op%Mk(iflavor,iflavor2,3)
1054   else
1055     op%Mk(iflavor,iflavor2,3) = 0.d0
1056     op%Mk(iflavor,iflavor2,2) = 0.d0
1057   endif
1058   op%setMk = op%setMk + 1
1059 !write(6,*) "mom1",op%Mk(iflavor,iflavor2,:)
1060 END SUBROUTINE GreenHyboffdiag_setMuD1

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setN [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_setN

FUNCTION

  impose number of electrons for this flavor

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green
  N=number of electrons

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

945 SUBROUTINE GreenHyboffdiag_setN(op,N)
946 
947 !Arguments ------------------------------------
948 
949 !This section has been created automatically by the script Abilint (TD).
950 !Do not modify the following lines by hand.
951 #undef ABI_FUNC
952 #define ABI_FUNC 'GreenHyboffdiag_setN'
953 !End of the abilint section
954 
955   TYPE(GreenHyboffdiag)    , INTENT(INOUT)    :: op
956   DOUBLE PRECISION  , INTENT(IN   )    :: N(op%nflavors)
957   INTEGER :: iflavor,iflavor2
958   COMPLEX(KIND=8) :: tmpoper
959 
960   IF ( op%set .EQV. .FALSE. ) &
961     CALL ERROR("GreenHyboffdiag_setN: green op%operator not set                ")
962   DO iflavor=1, op%nflavors
963    ! write(6,*) "iflavor",-N(iflavor)*op%signvaluemeas ,2*op%oper(op%samples,iflavor,iflavor),(N(iflavor)-1.d0)*op%signvaluemeas
964     ! the mulplication by signvaluemeas is necessary because N is
965     ! exactly the number of electrons in the flavor iflavor whereas
966     ! op%oper is not exact, because it still has to be divided by
967     ! signvaluemeas after the MPIREDUCE
968     op%oper(1,iflavor,iflavor) = (N(iflavor) - 1.d0)*op%signvaluemeas
969     op%oper(op%samples,iflavor,iflavor) = - N(iflavor)*op%signvaluemeas
970     !op%oper(op%samples,iflavor,iflavor) = 2*op%oper(op%samples,iflavor,iflavor)
971     !op%oper(1,iflavor,iflavor) = 2*op%oper(1,iflavor,iflavor)
972     DO iflavor2=1, op%nflavors
973       if(iflavor/=iflavor2) then
974               ! UNEXPLAINED but MANDATORY to have exact results for U=0 nspinor=4 with pawspnorb=0
975               ! Correction: The fact 2 is necessary for edge points because the points are at the
976               ! edges. 
977               ! It is of course necessary to fulfill exact results (U=0).
978     !tmpoper=(op%oper(op%samples,iflavor,iflavor2)-op%oper(1,iflavor,iflavor2))
979     !op%oper(op%samples,iflavor,iflavor2) = tmpoper
980     !op%oper(1,iflavor,iflavor2) = -tmpoper
981     !op%oper(op%samples,iflavor,iflavor2) = 2*op%oper(op%samples,iflavor,iflavor2)
982     !op%oper(1,iflavor,iflavor2) = 2*op%oper(1,iflavor,iflavor2)
983         op%oper(op%samples,iflavor,iflavor2) = 2*op%oper(op%samples,iflavor,iflavor2)
984         op%oper(1,iflavor,iflavor2) = 2*op%oper(1,iflavor,iflavor2)
985       endif
986     ENDDO
987   ENDDO
988 END SUBROUTINE GreenHyboffdiag_setN

ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setOperW [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyboffdiag_setOperW

FUNCTION

  set Green function in frequencies

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  op=Green
  Gomega=Input values

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

454 SUBROUTINE GreenHyboffdiag_setOperW(op, Gomega)
455 
456 !Arguments ------------------------------------
457 
458 !This section has been created automatically by the script Abilint (TD).
459 !Do not modify the following lines by hand.
460 #undef ABI_FUNC
461 #define ABI_FUNC 'GreenHyboffdiag_setOperW'
462 !End of the abilint section
463 
464   TYPE(GreenHyboffdiag)          , INTENT(INOUT) :: op
465   COMPLEX(KIND=8), DIMENSION(:,:,:), INTENT(IN   ) :: Gomega
466 !Loval variables ------------------------------
467   INTEGER :: tail
468 
469   tail = SIZE(Gomega,1)
470   IF ( .NOT. op%set ) &
471     CALL ERROR("GreenHyboffdiag_setOperW : Uninitialized GreenHyboffdiag structure")
472   IF ( ALLOCATED(op%oper_w) ) THEN
473     IF ( SIZE(op%oper_w) .NE. tail ) THEN
474       FREE(op%oper_w)
475       MALLOC(op%oper_w,(1:tail,op%nflavors,op%nflavors))
476     END IF
477   ELSE
478     MALLOC(op%oper_w,(1:tail,op%nflavors,op%nflavors))
479   END IF
480   op%oper_w(:,:,:) = Gomega(:,:,:)
481   op%Wmax = tail
482   op%setW = .TRUE.
483 END SUBROUTINE GreenHyboffdiag_setOperW

m_GreenHyboffdiag/GreenHyboffdiag [ Types ]

[ Top ] [ m_GreenHyboffdiag ] [ Types ]

NAME

  GreenHyboffdiag

FUNCTION

  This structured datatype contains the necessary data

COPYRIGHT

  Copyright (C) 2013 ABINIT group (J. Bieder)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

SOURCE

 76  TYPE GreenHyboffdiag
 77 
 78   LOGICAL :: set = .FALSE.
 79    ! True if variable of type GreenHyboffdiag is initialized
 80 
 81   LOGICAL :: setT = .FALSE.
 82    ! True if variable oper contains data
 83 
 84   LOGICAL :: setW = .FALSE.
 85    ! True if variable oper_w contains data
 86 
 87   LOGICAL :: have_MPI = .FALSE.
 88    ! True if MPI is used.
 89 
 90   INTEGER :: setMk = 0
 91    ! setMk=0 is moments for Fourier transform are not computed
 92 
 93   INTEGER :: samples
 94    ! samples=imaginary time slices (dmftqmc_l+1)
 95 
 96   INTEGER :: measurements    
 97    ! number of measurements for the Green's function
 98 
 99   INTEGER :: factor
100    ! if the move is not accepted, the statistic weight has to be
101    ! increased for the current configuration.
102 
103   INTEGER :: MY_COMM
104    ! MPI Communicator
105 
106   INTEGER :: size
107    ! size=1 
108 
109   INTEGER :: rank
110    ! rank=0 
111 
112   INTEGER :: Wmax
113    ! samples-1 if frequency Green's function
114 
115   INTEGER :: iTech
116    ! Precise if Frequency Green's function is computed or not
117 
118   INTEGER :: nflavors
119    ! Number of flavors 
120 
121   DOUBLE PRECISION :: beta
122    ! Inverse of temperature
123 
124   DOUBLE PRECISION :: inv_beta
125    ! Temperature
126 
127   DOUBLE PRECISION :: delta_t
128    ! 1/inv_dt
129 
130   DOUBLE PRECISION :: inv_dt
131    ! (samples-1)/beta
132   DOUBLE PRECISION :: signvaluemeas
133 
134   DOUBLE PRECISION :: signvalueold
135 
136   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: oper
137    ! oper(samples)
138 
139   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: omega
140    ! omega(Wmax) 
141 
142   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: Mk
143    ! Moments for FT
144 
145   COMPLEX(KIND=8)  , ALLOCATABLE, DIMENSION(:,:,:) :: oper_w 
146    ! Frequency Green's function
147 
148   COMPLEX(KIND=8)  , ALLOCATABLE, DIMENSION(:) :: oper_w_old
149    ! Old frequency Green's function (not used)
150 
151   TYPE(Vector)                            :: oper_old          
152    ! useless data
153 
154   TYPE(VectorInt)                         :: index_old          
155    ! useless data
156 
157   TYPE(MapHyb), ALLOCATABLE, DIMENSION(:,:)  :: map
158    ! value of time and Green's functions computed in GreenHyboffdiag_measHybrid
159    ! These values are used to fill op%oper in the same routine.
160 
161  END TYPE GreenHyboffdiag