TABLE OF CONTENTS
- ABINIT/m_GreenHyb
- ABINIT/m_GreenHyb/GreenHyb_backFourier
- ABINIT/m_GreenHyb/GreenHyb_clear
- ABINIT/m_GreenHyb/GreenHyb_destroy
- ABINIT/m_GreenHyb/GreenHyb_forFourier
- ABINIT/m_GreenHyb/GreenHyb_getHybrid
- ABINIT/m_GreenHyb/GreenHyb_init
- ABINIT/m_GreenHyb/GreenHyb_measHybrid
- ABINIT/m_GreenHyb/GreenHyb_print
- ABINIT/m_GreenHyb/GreenHyb_reset
- ABINIT/m_GreenHyb/GreenHyb_setMoments
- ABINIT/m_GreenHyb/GreenHyb_setMuD1
- ABINIT/m_GreenHyb/GreenHyb_setN
- ABINIT/m_GreenHyb/GreenHyb_setOperW
- m_GreenHyb/GreenHyb
ABINIT/m_GreenHyb [ Modules ]
NAME
m_GreenHyb
FUNCTION
Manage a green function for one orbital
COPYRIGHT
Copyright (C) 2013-2024 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 .
NOTES
SOURCE
22 #include "defs.h" 23 MODULE m_GreenHyb 24 USE m_Global 25 USE m_MatrixHyb 26 USE m_Vector 27 USE m_VectorInt 28 USE m_ListCdagC 29 USE m_MapHyb 30 #ifdef HAVE_MPI2 31 USE mpi 32 #endif 33 IMPLICIT NONE
ABINIT/m_GreenHyb/GreenHyb_backFourier [ Functions ]
NAME
GreenHyb_backFourier
FUNCTION
perform back fourier transform
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green dvgc=divergence parameter
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
688 SUBROUTINE GreenHyb_backFourier(this,dvgc) 689 690 691 #ifdef HAVE_MPI1 692 include 'mpif.h' 693 #endif 694 !Arguments ------------------------------------ 695 TYPE(GreenHyb) , INTENT(INOUT) :: this 696 DOUBLE PRECISION, OPTIONAL, INTENT(IN ) :: dvgc 697 !Local variables ------------------------------ 698 INTEGER :: itau 699 INTEGER :: iomega 700 INTEGER :: omegaSamples 701 INTEGER :: tauSamples 702 INTEGER :: tauBegin 703 INTEGER :: tauEnd 704 INTEGER :: delta 705 INTEGER :: residu 706 INTEGER, ALLOCATABLE, DIMENSION(:) :: counts 707 INTEGER, ALLOCATABLE, DIMENSION(:) :: displs 708 DOUBLE PRECISION :: A ! Correction factor 709 DOUBLE PRECISION :: inv_beta 710 DOUBLE PRECISION :: pi_invBeta 711 DOUBLE PRECISION :: two_invBeta 712 DOUBLE PRECISION :: minusDt 713 DOUBLE PRECISION :: minusOmegaTau 714 DOUBLE PRECISION :: omega 715 DOUBLE PRECISION :: minusTau 716 DOUBLE PRECISION :: sumTerm 717 DOUBLE PRECISION :: pi 718 DOUBLE PRECISION :: twoPi 719 DOUBLE PRECISION :: correction 720 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: Domega 721 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: A_omega 722 723 #if defined HAVE_MPI && !defined HAVE_MPI2_INPLACE 724 INTEGER :: my_count 725 DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:) :: oper_buf 726 #endif 727 728 IF ( this%set .EQV. .FALSE. ) & 729 CALL ERROR("GreenHyb_backFourier : Uninitialized GreenHyb structure") 730 IF ( this%setW .EQV. .FALSE. ) & 731 CALL ERROR("GreenHyb_backFourier : no G(iw)") 732 733 inv_beta = this%inv_beta 734 two_invBeta = 2.d0 * inv_beta 735 minusDt = - this%delta_t 736 omegaSamples = this%Wmax 737 tauSamples = this%samples-1 738 739 this%oper = 0.d0 740 741 pi = ACOS(-1.d0) 742 twoPi = 2.d0 * pi 743 pi_invBeta = pi * inv_beta 744 745 IF ( PRESENT(dvgc) ) THEN 746 A = dvgc 747 ELSE 748 A = AIMAG(this%oper_w(omegaSamples)) & ! A = \lim_\infty G(w) 749 *(2.d0*DBLE(omegaSamples)-1.d0) * pi_invBeta 750 END IF 751 752 correction = A*0.5d0 753 754 MALLOC(Domega,(1:omegaSamples)) 755 MALLOC(A_omega,(1:omegaSamples)) 756 Domega = (/ ((2.d0 * DBLE(iomega) - 1.d0)*pi_invbeta, iomega=1, omegaSamples) /) 757 A_omega = A / Domega 758 IF (this%have_MPI .EQV. .TRUE.) THEN 759 delta = tauSamples / this%size 760 residu = tauSamples - this%size*delta 761 IF ( this%rank .LT. this%size - residu ) THEN 762 tauBegin = 1 + this%rank*delta 763 tauEnd = (this%rank + 1)*delta 764 ELSE 765 ! tauBegin = (this%size-residu)*delta + 1 + (this%rank-this%size+residu)*(delta+1) 766 tauBegin = 1 + this%rank*(delta + 1) -this%size + residu 767 tauEnd = tauBegin + delta 768 END IF 769 MALLOC(counts,(1:this%size)) 770 MALLOC(displs,(1:this%size)) 771 counts = (/ (delta, iTau=1, this%size-residu), & 772 (delta+1, iTau=this%size-residu+1, this%size) /) 773 displs(1)=0 774 DO iTau = 2, this%size 775 displs(iTau) = displs(iTau-1) + counts (iTau-1) 776 END DO 777 ELSE 778 tauBegin = 1 779 tauEnd = tauSamples 780 END IF 781 DO iTau = tauBegin, tauEnd 782 minusTau = DBLE(itau -1) * minusDt 783 DO iomega = 1, omegaSamples 784 omega = Domega(iomega) 785 minusOmegaTau = MOD(omega*minusTau, TwoPi) 786 sumTerm = REAL(( this%oper_w(iomega) - CMPLX(0.d0, A_omega(iomega),8) ) & 787 * EXP( CMPLX(0.d0, minusOmegaTau, 8))) 788 this%oper(itau) = this%oper(itau) + sumTerm 789 END DO 790 this%oper(itau) = correction + two_invBeta*this%oper(itau) 791 END DO 792 IF ( this%have_MPI .EQV. .TRUE. ) THEN 793 ! rassembler les resultats 794 #ifdef HAVE_MPI 795 #if defined HAVE_MPI2_INPLACE 796 CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_PRECISION, & 797 this%oper, counts, displs, & 798 MPI_DOUBLE_PRECISION, this%MY_COMM, residu) 799 #else 800 my_count=tauBegin-tauEnd+1 801 MALLOC(oper_buf,(my_count)) 802 oper_buf(1:my_count)=this%oper(tauBegin:tauEnd) 803 CALL MPI_ALLGATHERV(oper_buf, my_count, MPI_DOUBLE_PRECISION, & 804 this%oper, counts, displs, & 805 MPI_DOUBLE_PRECISION, this%MY_COMM, residu) 806 FREE(oper_buf) 807 #endif 808 #endif 809 FREE(counts) 810 FREE(displs) 811 END IF 812 this%oper(tauSamples+1) = A - this%oper(1) !G(0+)-G(0-)=G(0+)+G(beta-)=A 813 this%setT = .TRUE. 814 FREE(Domega) 815 FREE(A_omega) 816 817 END SUBROUTINE GreenHyb_backFourier
ABINIT/m_GreenHyb/GreenHyb_clear [ Functions ]
NAME
GreenHyb_clear
FUNCTION
clear green function
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
237 SUBROUTINE GreenHyb_clear(this) 238 239 !Arguments ------------------------------------ 240 TYPE(GreenHyb) , INTENT(INOUT) :: this 241 242 !CALL Vector_clear(this%oper_old) 243 !CALL VectorInt_clear(this%index_old) 244 CALL MapHyb_clear(this%this) 245 this%measurements = 0 246 IF ( ALLOCATED(this%oper) ) & 247 this%oper = 0.d0 248 IF ( this%iTech .EQ. GREENHYB_OMEGA ) THEN 249 IF ( ALLOCATED(this%oper_w) ) & 250 this%oper_w = CMPLX(0.d0,0.d0,8) 251 IF ( ALLOCATED(this%oper_w_old) ) & 252 this%oper_w_old = CMPLX(0.d0,0.d0,8) 253 END IF 254 this%factor = 0 255 END SUBROUTINE GreenHyb_clear
ABINIT/m_GreenHyb/GreenHyb_destroy [ Functions ]
NAME
GreenHyb_destroy
FUNCTION
destroy green function
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1199 SUBROUTINE GreenHyb_destroy(this) 1200 1201 !Arguments ------------------------------------ 1202 TYPE(GreenHyb), INTENT(INOUT) :: this 1203 1204 this%set = .FALSE. 1205 this%setT = .FALSE. 1206 this%setW = .FALSE. 1207 this%samples = 0 1208 this%measurements = 0 1209 this%beta = 0.d0 1210 this%inv_beta = 0.d0 1211 this%inv_dt = 0.d0 1212 this%delta_t = 0.d0 1213 !CALL VectorInt_destroy(this%index_old) 1214 !CALL Vector_destroy(this%oper_old) 1215 CALL MapHyb_destroy(this%this) 1216 FREEIF(this%oper) 1217 FREEIF(this%oper_w) 1218 FREEIF(this%oper_w_old) 1219 FREEIF(this%omega) 1220 END SUBROUTINE GreenHyb_destroy
ABINIT/m_GreenHyb/GreenHyb_forFourier [ Functions ]
NAME
GreenHyb_forFourier
FUNCTION
perform forward fourier transform
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green Wmax=linear maximum frequency
OUTPUT
Gomega=Results for omega frequencies omega=ask frequencies
SIDE EFFECTS
NOTES
SOURCE
847 SUBROUTINE GreenHyb_forFourier(this, Gomega, omega, Wmax) 848 !Arguments ------------------------------------ 849 850 #ifdef HAVE_MPI1 851 include 'mpif.h' 852 #endif 853 TYPE(GreenHyb) , INTENT(INOUT) :: this 854 COMPLEX(KIND=8), DIMENSION(:), OPTIONAL, INTENT(INOUT) :: Gomega ! INOUT for MPI 855 COMPLEX(KIND=8), DIMENSION(:), OPTIONAL, INTENT(IN ) :: omega 856 INTEGER , OPTIONAL, INTENT(IN ) :: Wmax 857 INTEGER :: i 858 INTEGER :: j 859 INTEGER :: L 860 INTEGER :: Lspline 861 INTEGER :: Nom 862 INTEGER :: omegaBegin 863 INTEGER :: omegaEnd 864 INTEGER :: deltaw 865 INTEGER :: residu 866 INTEGER, ALLOCATABLE, DIMENSION(:) :: counts 867 INTEGER, ALLOCATABLE, DIMENSION(:) :: displs 868 DOUBLE PRECISION :: beta 869 DOUBLE PRECISION :: tau 870 DOUBLE PRECISION :: delta 871 DOUBLE PRECISION :: deltabis 872 DOUBLE PRECISION :: inv_delta 873 DOUBLE PRECISION :: inv_delta2 874 DOUBLE PRECISION :: omdeltabis 875 DOUBLE PRECISION :: tmp 876 DOUBLE PRECISION :: xpi 877 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: diag 878 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: diagL 879 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: lastR 880 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: lastC 881 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: XM 882 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: X2 883 DOUBLE PRECISION :: iw 884 COMPLEX(KIND=8) :: iwtau 885 COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:) :: Gwtmp 886 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: omegatmp 887 888 #if defined HAVE_MPI && !defined HAVE_MPI2_INPLACE 889 INTEGER :: my_count 890 COMPLEX(KIND=8), ALLOCATABLE , DIMENSION(:) :: Gwtmp_buf 891 #endif 892 893 IF ( this%set .EQV. .FALSE. ) & 894 CALL ERROR("GreenHyb_forFourier : Uninitialized GreenHyb structure") 895 IF ( this%setT .EQV. .FALSE. ) & 896 CALL ERROR("GreenHyb_forFourier : no G(tau)") 897 IF ( this%setMk .NE. 2 ) & 898 CALL WARNALL("GreenHyb_forFourier : green does not have moments ") 899 900 L = this%samples 901 902 xpi=acos(-1.d0) !!! XPI=PI 903 beta = this%beta 904 Nom = this%Wmax 905 IF ( PRESENT(Gomega) ) THEN 906 Nom = SIZE(Gomega) 907 !IF ( this%rank .EQ. 0 ) & 908 !write(6,*) "size Gomega", Nom 909 END IF 910 IF ( PRESENT(omega) ) THEN 911 IF ( PRESENT(Gomega) .AND. SIZE(omega) .NE. Nom ) THEN 912 CALL ERROR("GreenHyb_forFourier : sizes mismatch ") 913 !ELSE 914 !Nom = SIZE(omega) 915 END IF 916 END IF 917 IF ( .NOT. PRESENT(Gomega) .AND. .NOT. PRESENT(omega) ) THEN 918 IF ( PRESENT(Wmax) ) THEN 919 Nom=Wmax 920 ELSE 921 CALL ERROR("GreenHyb_forFourier : Missing argument Wmax") 922 END IF 923 END IF 924 925 IF ( ALLOCATED(this%oper_w) ) THEN 926 IF ( SIZE(this%oper_w) .NE. Nom ) THEN 927 FREE(this%oper_w) 928 MALLOC(this%oper_w,(1:Nom)) 929 END IF 930 ELSE 931 MALLOC(this%oper_w,(1:Nom)) 932 END IF 933 934 !write(6,*) "PRESENT(GOMEGA)", PRESENT(GOMEGA) 935 !write(6,*) "PRESENT(OMEGA)", PRESENT(OMEGA) 936 !call flush(6) 937 938 delta=this%delta_t 939 inv_delta = this%inv_dt 940 inv_delta2 = inv_delta*inv_delta 941 942 MALLOC(diagL,(L-1)) 943 MALLOC(lastR,(L-1)) 944 MALLOC(diag,(L)) 945 MALLOC(lastC,(L-1)) 946 947 !(cf Stoer) for the spline interpolation : 948 ! second derivatives XM solution of A*XM=B. 949 !A=(2.4.2.11) of Stoer&Bulirsch + 2 limit conditions 950 !The LU decomposition of A is known explicitly; 951 952 diag (1) = 4.d0 ! 1.d0 *4.d0 factor 4 added for conditionning 953 diagL(1) = 0.25d0 !1.d0/4.d0 954 lastR(1) = -0.5d0 ! -2.d0/4.d0 955 lastC(1) = 4.d0 ! 1.d0*4.d0 956 diag (2) = 4.d0 957 diagL(2) = 0.25d0 958 lastR(2) = -0.25d0 959 lastC(2) = -1.d0 960 961 DO i = 3, L-2 962 tmp = 4.d0 - diagL(i-1) 963 diagL(i) = 1.d0 / tmp 964 END DO 965 DO i = 3, L-2 966 diag (i) = 1.d0 / diagL(i) 967 lastR(i) = -(lastR(i-1)*diagL(i)) 968 lastC(i) = -(lastC(i-1)*diagL(i-1)) 969 END DO 970 971 tmp = 1.d0/diag(L-2) 972 diag (L-1) = 4.d0 - tmp 973 lastR(L-1) = (1.d0 - lastR(L-2))/ diag(L-1) 974 !diagL(L-1) = lastR(L-1) 975 diagL(L-1) = 0.d0 ! for the Lq=B resolution 976 !lastC(L-1) = 1.d0 - lastC(L-2)*diagL(L-1) ! equivalent to the next line 977 lastC(L-1) = 1.d0 - (lastC(L-2)*lastR(L-1)) ! True value 978 diag (L ) = 2.d0! - DOT_PRODUCT( lastR , lastC ) 979 tmp = 0.d0 980 DO i = 1, L-1 981 tmp = tmp + lastR(i)*lastC(i) 982 END DO 983 diag (L ) = diag (L ) - tmp 984 lastC(L-1) = lastC(L-1)-1.d0 ! 1 is removed for the u.XM=q resolution 985 986 ! construct the B vector from A.Xm=B 987 MALLOC(XM,(L)) 988 XM(1) = 4.d0*this%Mk(3) 989 XM(L) = (6.d0 * inv_delta) * ( this%Mk(2) - ( & 990 (this%oper(2)-this%oper(1)) + & 991 (this%oper(L)-this%oper(L-1)) ) * inv_delta ) 992 DO i = 2, L-1 993 XM(i) = (6.d0 * inv_delta2) * ( (this%oper(i+1) & 994 - 2.d0 * this%oper(i)) & 995 + this%oper(i-1) ) 996 END DO 997 998 ! Find second derivatives XM: Solve the system 999 ! SOLVING Lq= XM 1000 ! q = XM 1001 do j=1,L-1 1002 XM(j+1)=XM(j+1)-(diagL(j)*XM(j)) 1003 XM(L) =XM(L) -(lastR(j)*XM(j)) 1004 end do 1005 FREE(diagL) 1006 FREE(lastR) 1007 1008 ! SOLVING U.XM=q 1009 ! XM = q 1010 do j=L-1,2,-1 1011 XM(j+1) = XM(j+1) / diag(j+1) 1012 XM(j)= (XM(j)-(XM(L)*lastC(j)))-XM(j+1) 1013 end do 1014 XM(2) = XM(2) / diag(2) 1015 XM(1) = (XM(1)-XM(L)*lastC(1)) / diag(1) 1016 1017 FREE(diag) 1018 FREE(lastC) 1019 1020 Lspline = L-1 1021 MALLOC(X2,(1:Lspline+1)) ! We impose L = Nom 1022 !Construct L2 second derivative from known derivatives XM 1023 deltabis = beta / DBLE(Lspline) 1024 DO i = 1, Lspline 1025 tau = deltabis * DBLE(i-1) 1026 j = ((L-1)*(i-1))/Lspline + 1!INT(tau * inv_delta) + 1 1027 X2(i) = inv_delta * ( XM(j)*(DBLE(j)*delta - tau ) + XM(j+1)*(tau - DBLE(j-1)*delta) ) 1028 END DO 1029 X2(Lspline+1) = XM(L) 1030 FREE(XM) 1031 1032 IF ( this%have_MPI .EQV. .TRUE. ) THEN 1033 deltaw = Nom / this%size 1034 residu = Nom - this%size*deltaw 1035 IF ( this%rank .LT. this%size - residu ) THEN 1036 omegaBegin = 1 + this%rank*deltaw 1037 omegaEnd = (this%rank + 1)*deltaw 1038 ELSE 1039 ! tauBegin = (this%size-residu)*deltaw + 1 + (this%rank-this%size+residu)*(deltaw+1) 1040 omegaBegin = 1 + this%rank*(deltaw + 1) -this%size + residu 1041 omegaEnd = omegaBegin + deltaw 1042 END IF 1043 MALLOC(counts,(1:this%size)) 1044 MALLOC(displs,(1:this%size)) 1045 counts = (/ (deltaw, i=1, this%size-residu), & 1046 (deltaw+1, i=this%size-residu+1, this%size) /) 1047 displs(1)=0 1048 DO i = 2, this%size 1049 displs(i) = displs(i-1) + counts (i-1) 1050 END DO 1051 ELSE 1052 omegaBegin = 1 1053 omegaEnd = Nom 1054 END IF 1055 1056 this%Mk(1) = -1.d0 1057 MALLOC(omegatmp,(omegaBegin:omegaEnd)) 1058 IF ( PRESENT(omega) ) THEN 1059 omegatmp(omegaBegin:omegaEnd) = (/ (AIMAG(omega(i)),i=omegaBegin,omegaEnd) /) 1060 ELSE 1061 omegatmp(omegaBegin:omegaEnd) = (/ ((((2.d0*DBLE(i)-1.d0)*xpi)/Beta), i=omegaBegin,omegaEnd) /) 1062 END IF 1063 MALLOC(Gwtmp,(1:Nom)) 1064 DO i = omegaBegin, omegaEnd 1065 iw = omegatmp(i) 1066 omdeltabis = iw*deltabis 1067 Gwtmp(i)=CMPLX(0.d0,0.d0,8) 1068 DO j=2, Lspline ! We impose L+1 = Nom 1069 iwtau = CMPLX(0.d0,omdeltabis*DBLE(j-1),8) 1070 Gwtmp(i) = Gwtmp(i) + EXP(iwtau) * CMPLX((X2(j+1) + X2(j-1))-2.d0*X2(j),0.d0,8) 1071 END DO 1072 Gwtmp(i) = Gwtmp(i)/CMPLX(((iw*iw)*(iw*iw)*deltabis),0.d0,8) & 1073 1074 + CMPLX( & 1075 ( ((X2(2)-X2(1))+(X2(Lspline+1)-X2(Lspline)))/((iw*iw)*deltabis) -this%Mk(2) ) & 1076 /(iw*iw) & 1077 , & 1078 (this%Mk(1)-this%Mk(3)/(iw*iw))/iw & 1079 , 8) 1080 !+ CMPLX( (X2(2)-X2(1))+(X2(Lspline+1)-X2(Lspline)), 0.d0, 8 ) ) & 1081 ! / (((iw*iw)*(iw*iw))*CMPLX(deltabis,0.d0,8)) & 1082 !- CMPLX(this%Mk(1),0.d0,8)/iw & 1083 !+ CMPLX(this%Mk(2),0.d0,8)/(iw*iw) & 1084 !- CMPLX(this%Mk(3),0.d0,8)/((iw*iw)*iw) 1085 !IF ( this%rank .EQ. 0 ) write(12819,*) iw,gwtmp(i) 1086 END DO 1087 FREE(omegatmp) 1088 !call flush(12819) 1089 FREE(X2) 1090 IF ( this%have_MPI .EQV. .TRUE. ) THEN 1091 #ifdef HAVE_MPI 1092 #if defined HAVE_MPI2_INPLACE 1093 CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_COMPLEX, & 1094 Gwtmp , counts, displs, & 1095 MPI_DOUBLE_COMPLEX, this%MY_COMM, residu) 1096 #else 1097 my_count=omegaBegin-omegaEnd+1 1098 MALLOC(Gwtmp_buf,(my_count)) 1099 Gwtmp_buf(1:my_count)=Gwtmp(omegaBegin:omegaEnd) 1100 CALL MPI_ALLGATHERV(Gwtmp_buf, my_count, MPI_DOUBLE_COMPLEX, & 1101 Gwtmp , counts, displs, & 1102 MPI_DOUBLE_COMPLEX, this%MY_COMM, residu) 1103 FREE(Gwtmp_buf) 1104 #endif 1105 #endif 1106 FREE(counts) 1107 FREE(displs) 1108 END IF 1109 IF ( PRESENT(Gomega) ) THEN 1110 Gomega = Gwtmp 1111 END IF 1112 this%oper_w = Gwtmp 1113 this%setW = .TRUE. 1114 FREE(Gwtmp) 1115 END SUBROUTINE GreenHyb_forFourier
ABINIT/m_GreenHyb/GreenHyb_getHybrid [ Functions ]
NAME
GreenHyb_getHybrid
FUNCTION
reduce green function
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
521 SUBROUTINE GreenHyb_getHybrid(this) 522 523 !Arguments ------------------------------------ 524 TYPE(GreenHyb), INTENT(INOUT) :: this 525 526 IF ( this%set .EQV. .FALSE. ) & 527 CALL ERROR("GreenHyb_getHybrid : green operator not set ") 528 529 SELECT CASE(this%iTech) 530 CASE (GREENHYB_TAU) 531 this%oper = -(this%oper * this%inv_beta) / (DBLE(this%measurements) * this%delta_t) 532 this%setT = .TRUE. 533 CASE (GREENHYB_OMEGA) 534 this%oper_w = -(this%oper_w * this%inv_beta) / (DBLE(this%measurements) * this%delta_t) 535 this%setW = .TRUE. 536 CALL GreenHyb_backFourier(this) 537 END SELECT 538 539 END SUBROUTINE GreenHyb_getHybrid
ABINIT/m_GreenHyb/GreenHyb_init [ Functions ]
NAME
GreenHyb_init
FUNCTION
Initialize and allocate
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green samples=imaginary time slices beta=inverse temperature iTech=SHOULD NOT BE USED => BUGGY MY_COMM=mpi_communicator
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
126 SUBROUTINE GreenHyb_init(this, samples, beta,iTech,MY_COMM) 127 128 129 #ifdef HAVE_MPI1 130 include 'mpif.h' 131 #endif 132 !Arguments ------------------------------------ 133 TYPE(GreenHyb) , INTENT(INOUT) :: this 134 INTEGER , INTENT(IN ) :: samples 135 DOUBLE PRECISION, INTENT(IN ) :: beta 136 !INTEGER , INTENT(IN ) :: Wmax 137 INTEGER, OPTIONAL, INTENT(IN ) :: iTech 138 INTEGER, OPTIONAL, INTENT(IN ) :: MY_COMM 139 !Local variables ------------------------------ 140 INTEGER :: sp1 141 DOUBLE PRECISION :: dt 142 #ifdef HAVE_MPI 143 INTEGER :: ierr 144 #endif 145 146 IF ( PRESENT(MY_COMM)) THEN 147 #ifdef HAVE_MPI 148 this%have_MPI = .TRUE. 149 this%MY_COMM = MY_COMM 150 CALL MPI_Comm_rank(this%MY_COMM, this%rank, ierr) 151 CALL MPI_Comm_size(this%MY_COMM, this%size, ierr) 152 #else 153 CALL WARN("GreenHyb_init : MPI is not used ") 154 this%have_MPI = .FALSE. 155 this%MY_COMM = -1 156 this%rank = 0 157 this%size = 1 158 #endif 159 ELSE 160 this%have_MPI = .FALSE. 161 this%MY_COMM = -1 162 this%rank = 0 163 this%size = 1 164 END IF 165 166 sp1 = samples + 1 167 this%samples = sp1 168 this%measurements = 0 169 this%beta = beta 170 this%inv_beta = 1.d0 / beta 171 this%inv_dt = DBLE(samples) * this%inv_beta 172 dt = 1.d0 / this%inv_dt 173 this%delta_t = dt 174 !this%Wmax = Wmax 175 this%Wmax = -1 176 FREEIF(this%oper) 177 MALLOC(this%oper,(1:sp1)) 178 ! If we want to measure in frequences 179 ! let assume we first have "samples" frequences 180 IF ( PRESENT(iTech) ) THEN 181 this%iTech = iTech 182 SELECT CASE (this%iTech) 183 CASE (GREENHYB_TAU) ! omega 184 this%iTech = GREENHYB_TAU 185 CASE (GREENHYB_OMEGA) ! omega 186 this%Wmax = samples 187 FREEIF(this%oper_w) 188 MALLOC(this%oper_w,(1:this%Wmax)) 189 FREEIF(this%oper_w_old) 190 MALLOC(this%oper_w_old,(1:this%Wmax)) 191 this%oper_w = CMPLX(0.d0,0.d0,8) 192 this%oper_w_old = CMPLX(0.d0,0.d0,8) 193 FREEIF(this%omega) 194 MALLOC(this%omega,(1:this%Wmax)) 195 this%omega = (/ ((2.d0 * DBLE(sp1) - 1.d0)*ACOS(-1.d0)*this%inv_beta, sp1=1, this%Wmax) /) 196 END SELECT 197 ELSE 198 this%iTech = GREENHYB_TAU 199 END IF 200 ! end if 201 !CALL Vector_init(this%oper_old,10000) 202 !CALL VectorInt_init(this%index_old,10000) 203 CALL MapHyb_init(this%this,10000) 204 205 this%oper = 0.d0 206 this%set = .TRUE. 207 this%factor = 1 208 this%setMk = 0 209 this%Mk = 0.d0 210 END SUBROUTINE GreenHyb_init
ABINIT/m_GreenHyb/GreenHyb_measHybrid [ Functions ]
NAME
GreenHyb_measHybrid
FUNCTION
Measure Green function
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green Mthis=M this for the current flavor ListCdagC_1=list of all creator and annhilator operators updated=should we accumulate or not
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
372 SUBROUTINE GreenHyb_measHybrid(this, Mthis, ListCdagC_1, updated) 373 374 !Arguments ------------------------------------ 375 TYPE(GreenHyb) , INTENT(INOUT) :: this 376 TYPE(MatrixHyb) , INTENT(IN ) :: Mthis 377 TYPE(ListcdagC), INTENT(IN ) :: ListCdagC_1 378 LOGICAL , INTENT(IN ) :: updated 379 !Local variables ------------------------------ 380 INTEGER :: iC 381 INTEGER :: iCdag 382 INTEGER :: tail 383 !INTEGER :: index 384 INTEGER :: idx_old 385 INTEGER :: old_size 386 INTEGER :: omegaSamples 387 INTEGER :: iomega 388 DOUBLE PRECISION :: pi_invBeta 389 DOUBLE PRECISION :: mbeta_two 390 DOUBLE PRECISION :: beta 391 DOUBLE PRECISION :: beta_tc 392 DOUBLE PRECISION :: tcbeta_tc 393 DOUBLE PRECISION :: inv_dt 394 DOUBLE PRECISION :: tC 395 DOUBLE PRECISION :: tCdag 396 DOUBLE PRECISION :: time 397 DOUBLE PRECISION :: signe 398 DOUBLE PRECISION :: argument 399 !DOUBLE PRECISION :: taupi_invbeta 400 !COMPLEX(KIND=8) :: cargument 401 !COMPLEX(2*8) :: base_exp 402 !COMPLEX(2*8) :: increm_exp 403 404 IF ( this%set .EQV. .FALSE. ) & 405 CALL ERROR("GreenHyb_measHybrid : green operator not set ") 406 407 tail = ListCdagC_1%tail 408 IF ( tail .NE. Mthis%tail ) & 409 CALL ERROR("GreenHyb_measHybrid : ListCdagC & M unconsistent ") 410 411 IF ( updated .EQV. .TRUE. ) THEN ! NEW change in the configuration 412 ! FIXME SHOULD be much more faster 413 414 old_size = this%this%tail 415 SELECT CASE(this%iTech) 416 CASE (GREENHYB_TAU) 417 argument = DBLE(this%factor) 418 DO iC = 1, old_size 419 this%oper(this%this%listINT(iC)) = this%oper(this%this%listINT(iC)) + this%this%listDBLE(iC) * argument 420 END DO 421 this%measurements = this%measurements + this%factor 422 423 CALL MapHyb_setSize(this%this,tail*tail) 424 this%factor = 1 425 idx_old = 0 426 beta = this%beta 427 mbeta_two = -(beta*0.5d0) 428 inv_dt = this%inv_dt 429 ! WARNING time is not the time but just a temporary variable. 430 ! Index Time has been calculated previously and is in mat_tau 431 DO iC = 1, tail 432 tC = ListCdagC_1%list(iC,C_) 433 beta_tc = beta - tC 434 tcbeta_tc = tC * beta_tc 435 DO iCdag = 1, tail 436 tCdag = ListCdagC_1%list(iCdag,Cdag_) 437 time = tcbeta_tc - tCdag*beta_tc 438 439 !signe = SIGN(1.d0,time) 440 !time = time + (signe-1.d0)*mbeta_two 441 442 !signe = signe * SIGN(1.d0,beta-tC) 443 444 !signe = SIGN(1.d0,time) * SIGN(1.d0,beta-tC) 445 signe = SIGN(1.d0,time) 446 447 argument = signe*Mthis%mat(iCdag,iC) 448 449 !index = INT( ( time * inv_dt ) + 1.5d0 ) 450 !IF (index .NE. Mthis%mat_tau(iCdag,iC)) THEN 451 ! WRITE(*,*) index, Mthis%mat_tau(iCdag,iC) 452 !! CALL ERROR("Plantage") 453 !END IF 454 455 idx_old = idx_old + 1 456 this%this%listDBLE(idx_old) = argument 457 !this%this%listINT(idx_old) = index 458 this%this%listINT(idx_old) = Mthis%mat_tau(iCdag,iC) 459 END DO 460 END DO 461 CASE (GREENHYB_OMEGA) 462 omegaSamples = this%Wmax 463 argument = DBLE(this%factor) 464 DO iomega = 1, omegaSamples 465 this%oper_w(iomega) = this%oper_w(iomega) + this%oper_w_old(iomega) * argument 466 END DO 467 this%measurements = this%measurements + this%factor 468 469 this%factor = 1 470 beta = this%beta 471 mbeta_two = -(beta*0.5d0) 472 pi_invBeta = ACOS(-1.d0)/beta 473 DO iC = 1, tail 474 tC = ListCdagC_1%list(iC,C_) 475 DO iCdag = 1, tail 476 tCdag = ListCdagC_1%list(iCdag,Cdag_) 477 time = tC - tCdag 478 479 signe = SIGN(1.d0,time) 480 time = time + (signe-1.d0)*mbeta_two 481 signe = signe * SIGN(1.d0,beta-tC) 482 argument = signe*Mthis%mat(iCdag,iC) 483 484 DO iomega = 1, omegaSamples 485 !this%oper_w_old(iomega) = Mthis%mat_tau(iCdag,iC)*CMPLX(0.d0,argument) 486 this%oper_w_old(iomega) = EXP(CMPLX(0.d0,this%omega(iomega)*time,8))*CMPLX(0.d0,argument,8) 487 END DO 488 END DO 489 END DO 490 END SELECT 491 ELSE 492 this%factor = this%factor + 1 493 END IF 494 END SUBROUTINE GreenHyb_measHybrid
ABINIT/m_GreenHyb/GreenHyb_print [ Functions ]
NAME
GreenHyb_print
FUNCTION
print Green function
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green ostream=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1143 SUBROUTINE GreenHyb_print(this, ostream) 1144 1145 !Arguments ------------------------------------ 1146 TYPE(GreenHyb), INTENT(IN) :: this 1147 INTEGER, OPTIONAL , INTENT(IN) :: ostream 1148 !Local variables ------------------------------ 1149 INTEGER :: ostream_val 1150 INTEGER :: i 1151 INTEGER :: samples 1152 1153 1154 IF ( this%set .EQV. .FALSE. ) & 1155 CALL ERROR("GreenHyb_print : green this%operator not set ") 1156 1157 IF ( PRESENT(ostream) ) THEN 1158 ostream_val = ostream 1159 ELSE 1160 ostream_val = 66 1161 OPEN(UNIT=ostream_val,FILE="Green.dat") 1162 END IF 1163 1164 samples = this%samples 1165 1166 DO i = 1, samples 1167 WRITE(ostream_val,*) DBLE(i-1)*this%delta_t, this%oper(i) 1168 END DO 1169 1170 IF ( .NOT. PRESENT(ostream) ) & 1171 CLOSE(ostream_val) 1172 END SUBROUTINE GreenHyb_print
ABINIT/m_GreenHyb/GreenHyb_reset [ Functions ]
NAME
GreenHyb_reset
FUNCTION
reset green function
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
282 SUBROUTINE GreenHyb_reset(this) 283 284 !Arguments ------------------------------------ 285 TYPE(GreenHyb) , INTENT(INOUT) :: this 286 287 CALL GreenHyb_clear(this) 288 this%setMk = 0 289 this%Mk = 0.d0 290 this%setT = .FALSE. 291 this%setW = .FALSE. 292 END SUBROUTINE GreenHyb_reset
ABINIT/m_GreenHyb/GreenHyb_setMoments [ Functions ]
NAME
GreenHyb_setMoments
FUNCTION
Compute full moments
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green u1=interaction energi like u2=idem order2
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
646 SUBROUTINE GreenHyb_setMoments(this,u1,u2) 647 648 !Arguments ------------------------------------ 649 TYPE(GreenHyb) , INTENT(INOUT) :: this 650 DOUBLE PRECISION, INTENT(IN ) :: u1 651 DOUBLE PRECISION, INTENT(IN ) :: u2 652 653 this%Mk(1) = -1.d0 654 this%Mk(3) = this%Mk(3) - 2.d0*(this%Mk(2)*u1) 655 this%Mk(2) = this%Mk(2) + u1 656 this%Mk(3) = this%Mk(3) - u2 657 658 this%setMk = this%setMk + 1 659 660 END SUBROUTINE GreenHyb_setMoments
ABINIT/m_GreenHyb/GreenHyb_setMuD1 [ Functions ]
NAME
GreenHyb_setMuD1
FUNCTION
Set first moments for G
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green mu=energy level (irrespectige with fermi level) d1=firt moment of hybridization function
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
607 SUBROUTINE GreenHyb_setMuD1(this,mu,d1) 608 609 !Arguments ------------------------------------ 610 TYPE(GreenHyb) , INTENT(INOUT) :: this 611 DOUBLE PRECISION, INTENT(IN ) :: mu 612 DOUBLE PRECISION, INTENT(IN ) :: d1 613 614 this%Mk(3) = -d1-(mu*mu) 615 this%Mk(2) = -mu 616 this%setMk = this%setMk + 1 617 END SUBROUTINE GreenHyb_setMuD1
ABINIT/m_GreenHyb/GreenHyb_setN [ Functions ]
NAME
GreenHyb_setN
FUNCTION
impose number of electrons for this flavor
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green N=number of electrons
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
567 SUBROUTINE GreenHyb_setN(this,N) 568 569 !Arguments ------------------------------------ 570 TYPE(GreenHyb) , INTENT(INOUT) :: this 571 DOUBLE PRECISION, INTENT(IN ) :: N 572 573 IF ( this%set .EQV. .FALSE. ) & 574 CALL ERROR("GreenHyb_setN: green this%operator not set ") 575 this%oper(1) = N - 1.d0 576 this%oper(this%samples) = - N 577 END SUBROUTINE GreenHyb_setN
ABINIT/m_GreenHyb/GreenHyb_setOperW [ Functions ]
NAME
GreenHyb_setOperW
FUNCTION
set Green function in frequencies
COPYRIGHT
Copyright (C) 2013-2024 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
this=Green Gomega=Input values
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
320 SUBROUTINE GreenHyb_setOperW(this, Gomega) 321 322 !Arguments ------------------------------------ 323 TYPE(GreenHyb) , INTENT(INOUT) :: this 324 COMPLEX(KIND=8), DIMENSION(:), INTENT(IN ) :: Gomega 325 !Loval variables ------------------------------ 326 INTEGER :: tail 327 328 tail = SIZE(Gomega) 329 IF ( .NOT. this%set ) & 330 CALL ERROR("GreenHyb_setOperW : Uninitialized GreenHyb structure") 331 IF ( ALLOCATED(this%oper_w) ) THEN 332 IF ( SIZE(this%oper_w) .NE. tail ) THEN 333 FREE(this%oper_w) 334 MALLOC(this%oper_w,(1:tail)) 335 END IF 336 ELSE 337 MALLOC(this%oper_w,(1:tail)) 338 END IF 339 this%oper_w(:) = Gomega(:) 340 this%Wmax = tail 341 this%setW = .TRUE. 342 END SUBROUTINE GreenHyb_setOperW
m_GreenHyb/GreenHyb [ Types ]
[ Top ] [ m_GreenHyb ] [ Types ]
NAME
GreenHyb
FUNCTION
This structured datatype contains the necessary data
COPYRIGHT
Copyright (C) 2013-2024 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
54 TYPE, PUBLIC :: GreenHyb 55 LOGICAL _PRIVATE :: set = .FALSE. 56 LOGICAL _PRIVATE :: setT = .FALSE. 57 LOGICAL _PRIVATE :: setW = .FALSE. 58 LOGICAL _PRIVATE :: have_MPI = .FALSE. 59 INTEGER _PRIVATE :: setMk = 0 60 INTEGER _PRIVATE :: samples 61 INTEGER _PRIVATE :: measurements 62 INTEGER :: factor 63 INTEGER _PRIVATE :: MY_COMM 64 INTEGER _PRIVATE :: size 65 INTEGER _PRIVATE :: rank 66 INTEGER _PRIVATE :: Wmax 67 INTEGER _PRIVATE :: iTech 68 DOUBLE PRECISION _PRIVATE :: beta 69 DOUBLE PRECISION _PRIVATE :: inv_beta 70 DOUBLE PRECISION _PRIVATE :: delta_t 71 DOUBLE PRECISION _PRIVATE :: inv_dt 72 DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:) :: oper 73 DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:) _PRIVATE :: omega 74 DOUBLE PRECISION , DIMENSION(1:3) _PRIVATE :: Mk 75 COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:) _PRIVATE :: oper_w 76 COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:) _PRIVATE :: oper_w_old 77 TYPE(MapHyb) :: this 78 END TYPE GreenHyb