TABLE OF CONTENTS
- ABINIT/m_GreenHyboffdiag
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_backFourier
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_clear
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_destroy
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_forFourier
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_getHybrid
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_init
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_measHybrid
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_print
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_reset
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setMoments
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setMuD1
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setN
- ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setOperW
- m_GreenHyboffdiag/GreenHyboffdiag
ABINIT/m_GreenHyboffdiag [ Modules ]
NAME
m_GreenHyboffdiag
FUNCTION
Manage a green function for one orbital
COPYRIGHT
Copyright (C) 2013-2024 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
SOURCE
22 #include "defs.h" 23 MODULE m_GreenHyboffdiag 24 25 USE m_global 26 USE m_MatrixHyb 27 USE m_Vector 28 USE m_VectorInt 29 USE m_ListCdagC 30 USE m_MapHyb 31 #ifdef HAVE_MPI2 32 USE mpi 33 #endif 34 35 IMPLICIT NONE 36 37 public :: GreenHyboffdiag_init 38 public :: GreenHyboffdiag_reset 39 public :: GreenHyboffdiag_clear 40 public :: GreenHyboffdiag_setOperW 41 public :: GreenHyboffdiag_measHybrid 42 public :: GreenHyboffdiag_getHybrid 43 public :: GreenHyboffdiag_setN 44 public :: GreenHyboffdiag_setMuD1 45 public :: GreenHyboffdiag_setMoments 46 public :: GreenHyboffdiag_backFourier 47 public :: GreenHyboffdiag_forFourier 48 public :: GreenHyboffdiag_print 49 public :: GreenHyboffdiag_destroy 50 public :: nfourier3
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_backFourier [ Functions ]
NAME
GreenHyboffdiag_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
op=Green dvgc=divergence parameter
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1046 SUBROUTINE GreenHyboffdiag_backFourier(op,dvgc,func,hybri_limit,opt_hybri_limit) 1047 1048 use m_fstrings, only : int2char4 1049 1050 #ifdef HAVE_MPI1 1051 include 'mpif.h' 1052 #endif 1053 !Arguments ------------------------------------ 1054 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 1055 DOUBLE PRECISION, OPTIONAL, INTENT(IN ) :: dvgc 1056 CHARACTER(len=5) ,OPTIONAL, INTENT(IN) :: func 1057 COMPLEX(KIND=8), DIMENSION(op%nflavors,op%nflavors), OPTIONAL, INTENT(IN) :: hybri_limit 1058 INTEGER, OPTIONAL, INTENT(IN) :: opt_hybri_limit 1059 !Local variables ------------------------------ 1060 INTEGER :: itau 1061 INTEGER :: iomega 1062 INTEGER :: omegaSamples 1063 INTEGER :: tauSamples 1064 INTEGER :: tauBegin 1065 INTEGER :: tauEnd 1066 INTEGER :: delta 1067 INTEGER :: residu 1068 INTEGER :: iflavor1 1069 INTEGER :: iflavor2,unitnb !,unitnb1 1070 INTEGER, ALLOCATABLE, DIMENSION(:) :: counts 1071 INTEGER, ALLOCATABLE, DIMENSION(:) :: displs 1072 DOUBLE PRECISION :: A,AA ! Correction factor 1073 COMPLEX(KIND=8) :: B,BB ! Correction factor 1074 COMPLEX(KIND=8) :: C !,CC ! Correction factor 1075 DOUBLE PRECISION :: inv_beta 1076 DOUBLE PRECISION :: pi_invBeta 1077 DOUBLE PRECISION :: two_invBeta 1078 DOUBLE PRECISION :: minusDt 1079 DOUBLE PRECISION :: minusOmegaTau 1080 DOUBLE PRECISION :: omegaa 1081 DOUBLE PRECISION :: minusTau 1082 DOUBLE PRECISION :: sumTerm 1083 DOUBLE PRECISION :: pi 1084 DOUBLE PRECISION :: twoPi 1085 DOUBLE PRECISION :: correction 1086 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: Domega 1087 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: A_omega 1088 COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:) :: C_omega 1089 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: opertau 1090 CHARACTER(len=5) :: funct 1091 character(len=4) :: tag_proc 1092 character(len=30) :: tmpfil 1093 1094 #if defined HAVE_MPI && !defined HAVE_MPI2_INPLACE 1095 INTEGER :: my_count 1096 DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:) :: opertau_buf 1097 #endif 1098 1099 IF ( op%set .EQV. .FALSE. ) & 1100 CALL ERROR("GreenHyboffdiag_backFourier : Uninitialized GreenHyboffdiag structure") 1101 IF ( op%setW .EQV. .FALSE. ) & 1102 CALL ERROR("GreenHyboffdiag_backFourier : no G(iw)") 1103 1104 funct="hybri" 1105 if(present(func)) funct=func 1106 !sui!write(6,*) funct 1107 inv_beta = op%inv_beta 1108 two_invBeta = 2.d0 * inv_beta 1109 minusDt = - op%delta_t 1110 omegaSamples = op%Wmax 1111 tauSamples = op%samples-1 1112 pi = ACOS(-1.d0) 1113 twoPi = 2.d0 * pi 1114 pi_invBeta = pi * inv_beta 1115 !sui!write(6,*) "omegaSamples",omegaSamples 1116 MALLOC(Domega,(1:omegaSamples)) 1117 MALLOC(A_omega,(1:omegaSamples)) 1118 MALLOC(C_omega,(1:omegaSamples)) 1119 IF ( op%rank .EQ. 0 ) THEN 1120 !DO iflavor1 = 1, op%nflavors 1121 ! DO iflavor2 = 1, op%nflavors 1122 ! write(22236,*) "#",iflavor1,iflavor2 1123 ! do iomega=1,op%Wmax 1124 ! write(22236,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2)) 1125 ! enddo 1126 ! write(22236,*) 1127 ! ENDDO 1128 !ENDDO 1129 ENDIF 1130 1131 op%oper = 0.d0 1132 1133 DO iflavor1 = 1, op%nflavors 1134 DO iflavor2 = 1, op%nflavors 1135 ! -- compute limit of function G*(i\omega_n) 1136 if(funct=="hybri") then 1137 IF ( PRESENT(dvgc) ) THEN 1138 A = dvgc 1139 ELSE 1140 A = AIMAG(op%oper_w(omegaSamples,iflavor1,iflavor2))&! A = \lim_\infty Imag(G(iwn))*(wn) 1141 *(2.d0*DBLE(omegaSamples)-1.d0) * pi_invBeta 1142 AA = AIMAG(op%oper_w(omegaSamples-10,iflavor1,iflavor2))&! A = \lim_\infty Imag(G(iwn))*(wn) 1143 *(2.d0*DBLE(omegaSamples-10)-1.d0) * pi_invBeta 1144 B = op%oper_w(omegaSamples,iflavor1,iflavor2)&! A = \lim_\infty (G(iwn))*(iwn) 1145 *(2.d0*DBLE(omegaSamples)-1.d0) *pi_invBeta*cmplx(0.d0,1.d0,kind=8) 1146 BB = op%oper_w(omegaSamples-10,iflavor1,iflavor2)&! A = \lim_\infty (G(iwn))*(iwn) 1147 *(2.d0*DBLE(omegaSamples-10)-1.d0) *pi_invBeta*cmplx(0.d0,1.d0,kind=8) 1148 !sui!write(6,*) "B=",iflavor1,iflavor2,B,BB 1149 END IF 1150 else if(iflavor1==iflavor2.and.funct=="green") then 1151 A = -1.d0 1152 else if(iflavor1/=iflavor2.and.funct=="green") then 1153 A = 0.d0 1154 endif 1155 !sui!write(6,*) "A=",iflavor1,iflavor2,A,AA 1156 C=cmplx(-A,0.d0,kind=8) 1157 if(present(hybri_limit)) then 1158 if(present(opt_hybri_limit)) then 1159 if(opt_hybri_limit==1) C= (hybri_limit(iflavor1,iflavor2)) 1160 !sui!write(6,*) "C= ",C 1161 endif 1162 endif 1163 1164 1165 1166 ! -- correction on G(tau=0) is thus 1167 correction = -C*0.5d0 1168 1169 ! -- built frequency mesh 1170 Domega = (/ ((2.d0 * DBLE(iomega) - 1.d0)*pi_invbeta, iomega=1, omegaSamples) /) 1171 1172 ! -- built asymptotic function 1173 !if(present(hybri_limit)) then 1174 C_omega = C / (Domega*cmplx(0.d0,1.d0,kind=8)) 1175 !else 1176 ! A_omega = A / Domega 1177 ! C_omega=cmplx(0.d0,A_omega,kind=8) 1178 !endif 1179 !write(6,*) "AC 1",A_omega(2),C_omega(2) 1180 1181 IF ( op%rank .EQ. 0 ) THEN 1182 ! write(236,*) "#",iflavor1,iflavor2 1183 ! write(237,*) "#",iflavor1,iflavor2 1184 ! write(2236,*) "#",iflavor1,iflavor2 1185 ! write(2237,*) "#",iflavor1,iflavor2 1186 ! write(238,*) "#",iflavor1,iflavor2 1187 ! do iomega=1,op%Wmax 1188 ! write(2236,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2)) 1189 ! write(2237,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,imag(op%oper_w(iomega,iflavor1,iflavor2)) 1190 ! write(236,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2)-C_omega(iomega)) 1191 ! write(237,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,imag(op%oper_w(iomega,iflavor1,iflavor2)-C_omega(iomega)) 1192 ! write(238,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,C_omega(iomega) 1193 ! enddo 1194 ! write(236,*) 1195 ! write(237,*) 1196 ! write(2236,*) 1197 ! write(2237,*) 1198 ! write(238,*) 1199 ENDIF 1200 IF ( op%rank .EQ. 1 ) THEN 1201 ! write(22360,*) "#",iflavor1,iflavor2 1202 ! do iomega=1,op%Wmax 1203 ! write(22360,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2)) 1204 ! enddo 1205 ! write(22360,*) 1206 ENDIF 1207 !open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted') 1208 !write(unitnb,*) "#",iflavor1,iflavor2 1209 !do iomega=1,op%Wmax 1210 ! write(unitnb,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2)) 1211 !enddo 1212 !write(unitnb,*) 1213 ! -- built time mesh 1214 IF (op%have_MPI .EQV. .TRUE.) THEN 1215 delta = tauSamples / op%size 1216 residu = tauSamples - op%size*delta 1217 IF ( op%rank .LT. op%size - residu ) THEN 1218 tauBegin = 1 + op%rank*delta 1219 tauEnd = (op%rank + 1)*delta 1220 ELSE 1221 ! tauBegin = (op%size-residu)*delta + 1 + (op%rank-op%size+residu)*(delta+1) 1222 tauBegin = 1 + op%rank*(delta + 1) -op%size + residu 1223 tauEnd = tauBegin + delta 1224 END IF 1225 MALLOC(counts,(1:op%size)) 1226 MALLOC(displs,(1:op%size)) 1227 counts = (/ (delta, iTau=1, op%size-residu), & 1228 (delta+1, iTau=op%size-residu+1, op%size) /) 1229 displs(1)=0 1230 DO iTau = 2, op%size 1231 displs(iTau) = displs(iTau-1) + counts (iTau-1) 1232 END DO 1233 ELSE 1234 tauBegin = 1 1235 tauEnd = tauSamples 1236 END IF 1237 MALLOC(opertau,(1:tauSamples+1)) 1238 do iomega=1,omegaSamples 1239 ! write(6,*) iomega, imag(op%oper_w(iomega,iflavor1,iflavor2)), A_omega(iomega) ,"#diff" 1240 enddo 1241 unitnb=70000+op%rank 1242 call int2char4(op%rank,tag_proc) 1243 tmpfil = 'counts'//tag_proc 1244 ! open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted') 1245 ! write(unitnb,*) "#",iflavor1,iflavor2 1246 ! do itau=1,op%size 1247 ! write(unitnb,*) itau,counts(itau),displs(itau) 1248 ! enddo 1249 ! write(unitnb,*) 1250 1251 unitnb=10000+op%rank 1252 call int2char4(op%rank,tag_proc) 1253 tmpfil = 'oper_w'//tag_proc 1254 ! open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted') 1255 ! write(unitnb,*) "#",iflavor1,iflavor2,C 1256 ! ! C_omega et oper_w differents Domega identique. Est ce du a des 1257 !! ! diago differentes pour chaque procs dans qmc_prep_ctqmc 1258 !! do iomega=1,op%Wmax 1259 !! write(unitnb,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta,real(op%oper_w(iomega,iflavor1,iflavor2)),C_omega(iomega),Domega(iomega) 1260 ! enddo 1261 ! write(unitnb,*) 1262 1263 ! unitnb=40000+op%rank 1264 ! unitnb1=50000+op%rank 1265 ! call int2char4(op%rank,tag_proc) 1266 ! tmpfil = 'tauend'//tag_proc 1267 ! open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted') 1268 ! tmpfil = 'taubegin'//tag_proc 1269 ! open (unit=unitnb1,file=trim(tmpfil),status='unknown',form='formatted') 1270 ! write(unitnb,*) "#",iflavor1,iflavor2 1271 ! write(unitnb1,*) "#",iflavor1,iflavor2 1272 1273 ! -- compute Fourier transformation 1274 opertau=0.d0 1275 DO itau = tauBegin, tauEnd 1276 !DO itau = max(tauBegin-1,1), tauEnd 1277 minusTau = DBLE(itau -1) * minusDt 1278 DO iomega = 1, omegaSamples 1279 omegaa = Domega(iomega) 1280 minusOmegaTau = MOD(omegaa*minusTau, TwoPi) 1281 sumTerm = REAL(( op%oper_w(iomega,iflavor1,iflavor2) & 1282 - C_omega(iomega) ) & 1283 !- CMPLX(0.d0, A_omega(iomega),8) ) & 1284 * EXP( CMPLX(0.d0, minusOmegaTau, 8))) 1285 opertau(itau) = opertau(itau) + sumTerm 1286 ! Domega et minusomegatau identique MAIS oper_w different 1287 !write(unitnb,*) iomega,Domega(iomega),real(C_omega(iomega)),imag(C_omega(iomega)) 1288 !if(itau==tauend) then 1289 ! write(unitnb,*) iomega, sumTerm,opertau(itau),minusOmegaTau,op%oper_w(iomega,iflavor1,iflavor2),Domega(iomega) 1290 !endif 1291 !if(itau==max(tauBegin-1,1)) then 1292 ! write(unitnb1,*) iomega, sumTerm,opertau(itau),minusOmegaTau,op%oper_w(iomega,iflavor1,iflavor2),Domega(iomega) 1293 !endif 1294 1295 END DO 1296 ! if(itau==tauEnd) write(unitnb,*) 1297 ! if(itau==max(tauBegin-1,1)) write(unitnb1,*) 1298 ! if(iflavor1==iflavor2) then 1299 opertau(itau) = correction + two_invBeta*opertau(itau) 1300 ! if(itau==tauend) then 1301 ! write(unitnb,*) "final", opertau(itau),correction 1302 ! endif 1303 ! if(itau==max(tauBegin-1,1)) then 1304 ! write(unitnb1,*) "final",opertau(itau),correction 1305 ! endif 1306 !write(66666,*) itau, opertau(itau),correction 1307 ! else 1308 ! opertau(itau) = & 1309 ! two_invBeta*opertau(itau) 1310 ! endif 1311 END DO 1312 !unitnb=20000+op%rank 1313 !call int2char4(op%rank,tag_proc) 1314 !tmpfil = 'opertau'//tag_proc 1315 !open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted') 1316 !write(unitnb,*) "#",iflavor1,iflavor2,tauBegin,tauEnd 1317 !!do itau=tauBegin, tauEnd 1318 !do itau=1,tauSamples 1319 ! write(unitnb,*) itau,opertau(itau) 1320 !enddo 1321 !write(unitnb,*) 1322 !opertau(tauBegin-1)=0.d0 1323 !opertau(tauEnd+1)=0.d0 1324 !write(66666,*) 1325 1326 ! -- Gather 1327 IF ( op%have_MPI .EQV. .TRUE. ) THEN 1328 ! rassembler les resultats 1329 #ifdef HAVE_MPI 1330 #if defined HAVE_MPI2_INPLACE 1331 CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_PRECISION, & 1332 opertau, counts, displs, & 1333 MPI_DOUBLE_PRECISION, op%MY_COMM, residu) 1334 #else 1335 my_count=tauBegin-tauEnd+1 1336 MALLOC(opertau_buf,(my_count)) 1337 opertau_buf(1:my_count)=opertau(tauBegin:tauEnd) 1338 CALL MPI_ALLGATHERV(opertau_buf, my_count, MPI_DOUBLE_PRECISION, & 1339 opertau, counts, displs, & 1340 MPI_DOUBLE_PRECISION, op%MY_COMM, residu) 1341 FREE(opertau_buf) 1342 #endif 1343 #endif 1344 FREE(counts) 1345 FREE(displs) 1346 END IF 1347 ! unitnb=30000+op%rank 1348 ! call int2char4(op%rank,tag_proc) 1349 ! tmpfil = 'opertau_MPI_'//tag_proc 1350 ! open (unit=unitnb,file=trim(tmpfil),status='unknown',form='formatted') 1351 ! write(unitnb,*) "#",iflavor1,iflavor2 1352 ! do itau=tauBegin, tauEnd 1353 ! write(unitnb,*) itau,opertau(itau) 1354 ! enddo 1355 ! write(unitnb,*) 1356 1357 ! -- Add correction for discontinuity. 1358 ! if(iflavor1==iflavor2) then 1359 !G(0+)-G(0-)=G(0+)+G(beta-)=A 1360 opertau(tauSamples+1) = -real(C) - opertau(1) 1361 !sui!write(6,*) "BackFourier",opertau(tauSamples+1),opertau(1),real(C) 1362 1363 op%setT = .TRUE. 1364 ! endif 1365 op%oper(:,iflavor1,iflavor2)=opertau(:) 1366 FREE(opertau) 1367 END DO ! iflavor2 1368 END DO ! iflavor1 1369 ! -- End loop over flavors. 1370 1371 FREE(Domega) 1372 FREE(A_omega) 1373 FREE(C_omega) 1374 close(236) 1375 close(237) 1376 1377 END SUBROUTINE GreenHyboffdiag_backFourier
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_clear [ Functions ]
NAME
GreenHyboffdiag_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
op=Green
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
349 SUBROUTINE GreenHyboffdiag_clear(op) 350 351 !Arguments ------------------------------------ 352 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 353 INTEGER :: iflavor,iflavorbis 354 355 !CALL Vector_clear(op%oper_old) 356 !CALL VectorInt_clear(op%index_old) 357 do iflavor=1,op%nflavors 358 do iflavorbis=1,op%nflavors 359 CALL MapHyb_clear(op%map(iflavor,iflavorbis)) 360 enddo 361 enddo 362 op%measurements = 0 363 IF ( ALLOCATED(op%oper) ) & 364 op%oper = 0.d0 365 op%signvaluemeas = 0.d0 366 op%signvalueold = 1.d0 367 IF ( op%iTech .EQ. GREENHYB_OMEGA ) THEN 368 IF ( ALLOCATED(op%oper_w) ) & 369 op%oper_w = CMPLX(0.d0,0.d0,8) 370 IF ( ALLOCATED(op%oper_w_old) ) & 371 op%oper_w_old = CMPLX(0.d0,0.d0,8) 372 END IF 373 op%factor = 0 374 END SUBROUTINE GreenHyboffdiag_clear
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_destroy [ Functions ]
NAME
GreenHyboffdiag_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
op=Green
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1798 SUBROUTINE GreenHyboffdiag_destroy(op) 1799 1800 !Arguments ------------------------------------ 1801 TYPE(GreenHyboffdiag), INTENT(INOUT) :: op 1802 INTEGER :: iflavor,iflavorbis 1803 1804 op%set = .FALSE. 1805 op%setT = .FALSE. 1806 op%setW = .FALSE. 1807 op%samples = 0 1808 op%measurements = 0 1809 op%beta = 0.d0 1810 op%inv_beta = 0.d0 1811 op%inv_dt = 0.d0 1812 op%delta_t = 0.d0 1813 CALL VectorInt_destroy(op%index_old) 1814 CALL Vector_destroy(op%oper_old) 1815 do iflavor=1,op%nflavors 1816 do iflavorbis=1,op%nflavors 1817 !sui!write(6,*) "test",iflavor,iflavorbis 1818 CALL MapHyb_destroy(op%map(iflavor,iflavorbis)) 1819 enddo 1820 enddo 1821 DT_FREEIF(op%map) 1822 FREEIF(op%oper) 1823 FREEIF(op%Mk) 1824 FREEIF(op%oper_w) 1825 FREEIF(op%oper_w_old) 1826 FREEIF(op%omega) 1827 END SUBROUTINE GreenHyboffdiag_destroy
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_forFourier [ Functions ]
NAME
GreenHyboffdiag_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
op=Green Wmax=linear maximum frequency
OUTPUT
Gomega=Results for omega frequencies omega=ask frequencies
SIDE EFFECTS
NOTES
SOURCE
1407 SUBROUTINE GreenHyboffdiag_forFourier(op, Gomega, omega, Wmax) 1408 !Arguments ------------------------------------ 1409 1410 #ifdef HAVE_MPI1 1411 include 'mpif.h' 1412 #endif 1413 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 1414 COMPLEX(KIND=8), DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT) :: Gomega ! INOUT for MPI 1415 COMPLEX(KIND=8), DIMENSION(:), OPTIONAL, INTENT(IN ) :: omega 1416 INTEGER , OPTIONAL, INTENT(IN ) :: Wmax 1417 INTEGER :: i 1418 INTEGER :: j 1419 INTEGER :: iflavor1 1420 INTEGER :: iflavor2 1421 INTEGER :: nflavors 1422 INTEGER :: L 1423 INTEGER :: Lspline 1424 INTEGER :: Nom 1425 INTEGER :: omegaBegin 1426 INTEGER :: omegaEnd 1427 INTEGER :: deltaw 1428 INTEGER :: residu 1429 INTEGER, ALLOCATABLE, DIMENSION(:) :: counts 1430 INTEGER, ALLOCATABLE, DIMENSION(:) :: displs 1431 DOUBLE PRECISION :: beta 1432 DOUBLE PRECISION :: tau 1433 DOUBLE PRECISION :: delta 1434 DOUBLE PRECISION :: deltabis 1435 DOUBLE PRECISION :: inv_delta 1436 DOUBLE PRECISION :: inv_delta2 1437 DOUBLE PRECISION :: omdeltabis 1438 DOUBLE PRECISION :: tmp 1439 DOUBLE PRECISION :: xpi 1440 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: diag 1441 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: diagL 1442 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: lastR 1443 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: lastC 1444 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: XM 1445 DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: X2 1446 DOUBLE PRECISION :: iw 1447 COMPLEX(KIND=8) :: iwtau 1448 COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:) :: Gwtmp 1449 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: omegatmp 1450 1451 #if defined HAVE_MPI && !defined HAVE_MPI2_INPLACE 1452 INTEGER :: my_count 1453 COMPLEX(KIND=8), ALLOCATABLE , DIMENSION(:) :: Gwtmp_buf 1454 #endif 1455 1456 nflavors=op%nflavors 1457 1458 !sui!write(6,*) " Fourier transformation begin" 1459 1460 IF ( op%set .EQV. .FALSE. ) & 1461 CALL ERROR("GreenHyboffdiag_forFourier : Uninitialized GreenHyboffdiag structure") 1462 IF ( op%setT .EQV. .FALSE. ) & 1463 CALL ERROR("GreenHyboffdiag_forFourier : no G(tau)") 1464 !write(6,*) "op%setMk=", op%setMk 1465 IF ( op%setMk .NE. 2*nflavors*nflavors ) & 1466 CALL WARNALL("GreenHyboffdiag_forFourier : green does not have moments ") 1467 1468 L = op%samples 1469 1470 xpi=acos(-1.d0) !!! XPI=PI 1471 beta = op%beta 1472 Nom = op%Wmax 1473 IF ( PRESENT(Gomega) ) THEN 1474 Nom = SIZE(Gomega,1) 1475 !IF ( op%rank .EQ. 0 ) & 1476 !!write(6,*) "size Gomega", Nom 1477 END IF 1478 IF ( PRESENT(omega) ) THEN 1479 IF ( PRESENT(Gomega) .AND. SIZE(omega) .NE. Nom ) THEN 1480 CALL ERROR("GreenHyboffdiag_forFourier : sizes mismatch ") 1481 !ELSE 1482 !Nom = SIZE(omega) 1483 END IF 1484 END IF 1485 IF ( .NOT. PRESENT(Gomega) .AND. .NOT. PRESENT(omega) ) THEN 1486 IF ( PRESENT(Wmax) ) THEN 1487 Nom=Wmax 1488 ELSE 1489 CALL ERROR("GreenHyboffdiag_forFourier : Missing argument Wmax") 1490 END IF 1491 END IF 1492 1493 !!IF ( ALLOCATED(op%oper_w) ) THEN 1494 !! IF ( SIZE(op%oper_w,1) .NE. Nom ) THEN 1495 !! FREE(op%oper_w) 1496 !! MALLOC(op%oper_w,(1:Nom,nflavors,nflavors)) 1497 !! END IF 1498 !!ELSE 1499 !! MALLOC(op%oper_w,(1:Nom,nflavors,nflavors)) 1500 !!END IF 1501 1502 !!write(6,*) "PRESENT(GOMEGA)", PRESENT(GOMEGA) 1503 !!write(6,*) "PRESENT(OMEGA)", PRESENT(OMEGA) 1504 !call flush(6) 1505 1506 delta=op%delta_t 1507 inv_delta = op%inv_dt 1508 inv_delta2 = inv_delta*inv_delta 1509 1510 MALLOC(diagL,(L-1)) 1511 MALLOC(lastR,(L-1)) 1512 MALLOC(diag,(L)) 1513 MALLOC(lastC,(L-1)) 1514 1515 !(cf Stoer) for the spline interpolation : 1516 ! second derivatives XM solution of A*XM=B. 1517 !A=(2.4.2.11) of Stoer&Bulirsch + 2 limit conditions 1518 !The LU decomposition of A is known explicitly; 1519 1520 diag (1) = 4.d0 ! 1.d0 *4.d0 factor 4 added for conditionning 1521 diagL(1) = 0.25d0 !1.d0/4.d0 1522 lastR(1) = -0.5d0 ! -2.d0/4.d0 1523 lastC(1) = 4.d0 ! 1.d0*4.d0 1524 diag (2) = 4.d0 1525 diagL(2) = 0.25d0 1526 lastR(2) = -0.25d0 1527 lastC(2) = -1.d0 1528 1529 DO i = 3, L-2 1530 tmp = 4.d0 - diagL(i-1) 1531 diagL(i) = 1.d0 / tmp 1532 END DO 1533 DO i = 3, L-2 1534 diag (i) = 1.d0 / diagL(i) 1535 lastR(i) = -(lastR(i-1)*diagL(i)) 1536 lastC(i) = -(lastC(i-1)*diagL(i-1)) 1537 END DO 1538 1539 tmp = 1.d0/diag(L-2) 1540 diag (L-1) = 4.d0 - tmp 1541 lastR(L-1) = (1.d0 - lastR(L-2))/ diag(L-1) 1542 !diagL(L-1) = lastR(L-1) 1543 diagL(L-1) = 0.d0 ! for the Lq=B resolution 1544 !lastC(L-1) = 1.d0 - lastC(L-2)*diagL(L-1) ! equivalent to the next line 1545 lastC(L-1) = 1.d0 - (lastC(L-2)*lastR(L-1)) ! True value 1546 diag (L ) = 2.d0! - DOT_PRODUCT( lastR , lastC ) 1547 tmp = 0.d0 1548 DO i = 1, L-1 1549 tmp = tmp + lastR(i)*lastC(i) 1550 END DO 1551 diag (L ) = diag (L ) - tmp 1552 lastC(L-1) = lastC(L-1)-1.d0 ! 1 is removed for the u.XM=q resolution 1553 1554 MALLOC(XM,(L)) 1555 MALLOC(Gwtmp,(1:Nom)) 1556 1557 Lspline = L-1 1558 MALLOC(X2,(1:Lspline+1)) ! We impose L = Nom 1559 1560 IF ( op%have_MPI .EQV. .TRUE. ) THEN 1561 deltaw = Nom / op%size 1562 residu = Nom - op%size*deltaw 1563 IF ( op%rank .LT. op%size - residu ) THEN 1564 omegaBegin = 1 + op%rank*deltaw 1565 omegaEnd = (op%rank + 1)*deltaw 1566 ELSE 1567 ! tauBegin = (op%size-residu)*deltaw + 1 + (op%rank-op%size+residu)*(deltaw+1) 1568 omegaBegin = 1 + op%rank*(deltaw + 1) -op%size + residu 1569 omegaEnd = omegaBegin + deltaw 1570 END IF 1571 MALLOC(counts,(1:op%size)) 1572 MALLOC(displs,(1:op%size)) 1573 counts = (/ (deltaw, i=1, op%size-residu), & 1574 (deltaw+1, i=op%size-residu+1, op%size) /) 1575 displs(1)=0 1576 DO i = 2, op%size 1577 displs(i) = displs(i-1) + counts (i-1) 1578 END DO 1579 ELSE 1580 omegaBegin = 1 1581 omegaEnd = Nom 1582 END IF 1583 1584 ! op%Mk(iflavor1,iflavor2,1) = 0.d0 1585 ! DO iflavor1 = 1, nflavors 1586 ! op%Mk(iflavor1,iflavor1,1) = -1.d0 1587 ! ENDDO 1588 ! op%Mk(:,:,3) = 0.d0 1589 1590 MALLOC(omegatmp,(omegaBegin:omegaEnd)) 1591 IF ( PRESENT(omega) ) THEN 1592 omegatmp(omegaBegin:omegaEnd) = (/ (AIMAG(omega(i)),i=omegaBegin,omegaEnd) /) 1593 ELSE 1594 omegatmp(omegaBegin:omegaEnd) = (/ ((((2.d0*DBLE(i)-1.d0)*xpi)/Beta), i=omegaBegin,omegaEnd) /) 1595 END IF 1596 1597 DO iflavor1 = 1, nflavors 1598 DO iflavor2 = 1, nflavors 1599 ! write(6,*) " Moments:",op%Mk(iflavor1,iflavor2,:),iflavor1,iflavor2 1600 1601 ! construct the B vector from A.Xm=B 1602 XM(1) = 4.d0*op%Mk(iflavor1,iflavor2,3) 1603 XM(L) = (6.d0 * inv_delta) * ( op%Mk(iflavor1,iflavor2,2) - ( & 1604 (op%oper(2,iflavor1,iflavor2)-op%oper(1,iflavor1,iflavor2)) + & 1605 (op%oper(L,iflavor1,iflavor2)-op%oper(L-1,iflavor1,iflavor2)) ) * inv_delta ) 1606 ! built generic second derivative of oper 1607 !sui!write(6,*) "XM 1 L",XM(1),XM(L),op%Mk(iflavor1,iflavor2,2),op%Mk(iflavor1,iflavor2,3) 1608 DO i = 2, L-1 1609 XM(i) = (6.d0 * inv_delta2) * ( (op%oper(i+1,iflavor1,iflavor2) & 1610 - 2.d0 * op%oper(i,iflavor1,iflavor2)) & 1611 + op%oper(i-1,iflavor1,iflavor2) ) 1612 !sui!write(6,*) "XM",i,XM(i),op%oper(i,iflavor1,iflavor2) 1613 END DO 1614 1615 ! Find second derivatives XM: Solve the system 1616 ! SOLVING Lq= XM 1617 ! q = XM 1618 do j=1,L-1 1619 XM(j+1)=XM(j+1)-(diagL(j)*XM(j)) 1620 XM(L) =XM(L) -(lastR(j)*XM(j)) 1621 end do 1622 1623 1624 ! SOLVING U.XM=q 1625 ! XM = q 1626 do j=L-1,2,-1 1627 XM(j+1) = XM(j+1) / diag(j+1) 1628 XM(j)= (XM(j)-(XM(L)*lastC(j)))-XM(j+1) 1629 end do 1630 XM(2) = XM(2) / diag(2) 1631 XM(1) = (XM(1)-XM(L)*lastC(1)) / diag(1) 1632 1633 1634 1635 !Construct L2 second derivative from known derivatives XM 1636 deltabis = beta / DBLE(Lspline) 1637 DO i = 1, Lspline 1638 tau = deltabis * DBLE(i-1) 1639 j = ((L-1)*(i-1))/Lspline + 1!INT(tau * inv_delta) + 1 1640 X2(i) = inv_delta * ( XM(j)*(DBLE(j)*delta - tau ) + XM(j+1)*(tau - DBLE(j-1)*delta) ) 1641 END DO 1642 X2(Lspline+1) = XM(L) 1643 1644 1645 DO i = omegaBegin, omegaEnd 1646 iw = omegatmp(i) 1647 omdeltabis = iw*deltabis 1648 Gwtmp(i)=CMPLX(0.d0,0.d0,8) 1649 DO j=2, Lspline ! We impose L+1 = Nom 1650 iwtau = CMPLX(0.d0,omdeltabis*DBLE(j-1),8) 1651 Gwtmp(i) = Gwtmp(i) + EXP(iwtau) * CMPLX((X2(j+1) + X2(j-1))-2.d0*X2(j),0.d0,8) 1652 !write(6,*) "ww",i,j,Gwtmp(i),X2(j),iwtau 1653 END DO 1654 Gwtmp(i) = Gwtmp(i)/CMPLX(((iw*iw)*(iw*iw)*deltabis),0.d0,8) & 1655 + CMPLX( ( ((X2(2)-X2(1))+(X2(Lspline+1)-X2(Lspline)))/((iw*iw)*deltabis) -op%Mk(iflavor1,iflavor2,2) ) & 1656 /(iw*iw) , (op%Mk(iflavor1,iflavor2,1)-op%Mk(iflavor1,iflavor2,3)/(iw*iw))/iw , 8) 1657 !+ CMPLX( (X2(2)-X2(1))+(X2(Lspline+1)-X2(Lspline)), 0.d0, 8 ) ) & 1658 ! / (((iw*iw)*(iw*iw))*CMPLX(deltabis,0.d0,8)) & 1659 !- CMPLX(op%Mk(1),0.d0,8)/iw & 1660 !+ CMPLX(op%Mk(2),0.d0,8)/(iw*iw) & 1661 !- CMPLX(op%Mk(3),0.d0,8)/((iw*iw)*iw) 1662 !IF ( op%rank .EQ. 0 ) write(12819,*) iw,gwtmp(i) 1663 END DO 1664 !call flush(12819) 1665 IF ( op%have_MPI .EQV. .TRUE. ) THEN 1666 #ifdef HAVE_MPI 1667 #if defined HAVE_MPI2_INPLACE 1668 CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_COMPLEX, & 1669 Gwtmp , counts, displs, & 1670 MPI_DOUBLE_COMPLEX, op%MY_COMM, residu) 1671 #else 1672 my_count=omegaBegin-omegaEnd+1 1673 MALLOC(Gwtmp_buf,(my_count)) 1674 Gwtmp_buf(1:my_count)=Gwtmp(omegaBegin:omegaEnd) 1675 CALL MPI_ALLGATHERV(Gwtmp_buf, my_count, MPI_DOUBLE_COMPLEX, & 1676 Gwtmp , counts, displs, & 1677 MPI_DOUBLE_COMPLEX, op%MY_COMM, residu) 1678 FREE(Gwtmp_buf) 1679 #endif 1680 #endif 1681 END IF 1682 IF ( PRESENT(Gomega) ) THEN 1683 Gomega(:,iflavor1,iflavor2) = Gwtmp(:) 1684 END IF 1685 op%setW = .TRUE. 1686 ENDDO ! iflavor1 1687 ENDDO ! iflavor2 1688 !!op%oper_w=Gomega 1689 do iflavor1=1,nflavors 1690 !sui!write(6,*) iflavor1 1691 do i=1,Nom 1692 !write(6,*) "w",i,op%oper_w(i,iflavor1,iflavor1) 1693 enddo 1694 enddo 1695 FREE(Gwtmp) 1696 FREE(diagL) 1697 FREE(lastR) 1698 FREE(diag) 1699 FREE(lastC) 1700 FREE(XM) 1701 FREE(omegatmp) 1702 FREE(X2) 1703 FREE(counts) 1704 FREE(displs) 1705 1706 END SUBROUTINE GreenHyboffdiag_forFourier
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_getHybrid [ Functions ]
NAME
GreenHyboffdiag_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
op=Green
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
805 SUBROUTINE GreenHyboffdiag_getHybrid(op) 806 807 !Arguments ------------------------------------ 808 TYPE(GreenHyboffdiag), INTENT(INOUT) :: op 809 810 IF ( op%set .EQV. .FALSE. ) & 811 CALL ERROR("GreenHyboffdiag_getHybrid : green operator not set ") 812 813 SELECT CASE(op%iTech) 814 CASE (GREENHYB_TAU) 815 op%oper = -(op%oper * op%inv_beta) / (DBLE(op%measurements) * op%delta_t) 816 !sui!write(6,*) "measurements",op%measurements,op%delta_t,op%inv_beta 817 !sui!write(6,*) "signevaluemeas meas",op%signvaluemeas,op%measurements 818 op%signvaluemeas = op%signvaluemeas / DBLE(op%measurements) 819 ! print*, "op%oper",op%oper(1,1,1) 820 !sui!write(6,*) "signevaluemeas/meas",op%signvaluemeas 821 ! print*, "signevaluemeas/meas",op%signvaluemeas 822 op%setT = .TRUE. 823 CASE (GREENHYB_OMEGA) 824 op%oper_w = -(op%oper_w * op%inv_beta) / (DBLE(op%measurements) * op%delta_t) 825 op%setW = .TRUE. 826 CALL GreenHyboffdiag_backFourier(op) 827 END SELECT 828 829 END SUBROUTINE GreenHyboffdiag_getHybrid
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_init [ Functions ]
NAME
GreenHyboffdiag_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
op=Green samples=imaginary time slices beta=inverse temperature iTech=SHOULD NOT BE USED => BUGGY MY_COMM=mpi_communicator
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
188 SUBROUTINE GreenHyboffdiag_init(op, samples, beta,nflavors,iTech,MY_COMM) 189 190 191 #ifdef HAVE_MPI1 192 include 'mpif.h' 193 #endif 194 !Arguments ------------------------------------ 195 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 196 INTEGER , INTENT(IN ) :: samples 197 DOUBLE PRECISION, INTENT(IN ) :: beta 198 INTEGER , INTENT(IN ) :: nflavors 199 !INTEGER , INTENT(IN ) :: Wmax 200 INTEGER, OPTIONAL, INTENT(IN ) :: iTech 201 INTEGER, OPTIONAL, INTENT(IN ) :: MY_COMM 202 !Local variables ------------------------------ 203 INTEGER :: iflavor,iflavorbis,sp1 204 DOUBLE PRECISION :: dt 205 #ifdef HAVE_MPI 206 INTEGER :: ierr 207 #endif 208 209 IF ( PRESENT(MY_COMM)) THEN 210 #ifdef HAVE_MPI 211 op%have_MPI = .TRUE. 212 op%MY_COMM = MY_COMM 213 CALL MPI_Comm_rank(op%MY_COMM, op%rank, ierr) 214 CALL MPI_Comm_size(op%MY_COMM, op%size, ierr) 215 #else 216 CALL WARN("GreenHyboffdiag_init : MPI is not used ") 217 op%have_MPI = .FALSE. 218 op%MY_COMM = -1 219 op%rank = 0 220 op%size = 1 221 #endif 222 ELSE 223 op%have_MPI = .FALSE. 224 op%MY_COMM = -1 225 op%rank = 0 226 op%size = 1 227 END IF 228 229 sp1 = samples + 1 230 op%samples = sp1 231 op%measurements = 0 232 op%nflavors = nflavors 233 op%beta = beta 234 op%inv_beta = 1.d0 / beta 235 op%inv_dt = DBLE(samples) * op%inv_beta 236 dt = 1.d0 / op%inv_dt 237 op%delta_t = dt 238 !op%Wmax = Wmax 239 op%Wmax = -1 240 FREEIF(op%Mk) 241 MALLOC(op%Mk,(nflavors,nflavors,3)) 242 FREEIF(op%oper) 243 MALLOC(op%oper,(sp1,nflavors,nflavors)) 244 ! If we want to measure in frequences 245 ! let assume we first have "samples" frequences 246 IF ( PRESENT(iTech) ) THEN 247 op%iTech = iTech 248 SELECT CASE (op%iTech) 249 CASE (GREENHYB_TAU) ! omega 250 op%iTech = GREENHYB_TAU 251 CASE (GREENHYB_OMEGA) ! omega 252 op%Wmax = samples 253 FREEIF(op%oper_w) 254 MALLOC(op%oper_w,(1:op%Wmax,nflavors,nflavors)) 255 FREEIF(op%oper_w_old) 256 MALLOC(op%oper_w_old,(1:op%Wmax)) 257 op%oper_w = CMPLX(0.d0,0.d0,8) 258 op%oper_w_old = CMPLX(0.d0,0.d0,8) 259 FREEIF(op%omega) 260 MALLOC(op%omega,(1:op%Wmax)) 261 op%omega = (/ ((2.d0 * DBLE(sp1) - 1.d0)*ACOS(-1.d0)*op%inv_beta, sp1=1, op%Wmax) /) 262 END SELECT 263 ELSE 264 op%iTech = GREENHYB_TAU 265 END IF 266 ! end if 267 CALL Vector_init(op%oper_old,10000) 268 CALL VectorInt_init(op%index_old,10000) 269 DT_FREEIF(op%map) 270 MALLOC(op%map,(nflavors,nflavors)) 271 do iflavor=1,nflavors 272 do iflavorbis=1,nflavors 273 CALL MapHyb_init(op%map(iflavor,iflavorbis),10000) 274 enddo 275 enddo 276 277 op%oper = 0.d0 278 op%signvaluemeas = 0.d0 279 op%signvalueold = 0.d0 280 op%set = .TRUE. 281 op%factor = 1 282 op%setMk = 0 283 op%Mk = 0.d0 284 END SUBROUTINE GreenHyboffdiag_init
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_measHybrid [ Functions ]
NAME
GreenHyboffdiag_measHybrid
FUNCTION
Measure Green's 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
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
SOURCE
454 SUBROUTINE GreenHyboffdiag_measHybrid(op, Mmatrix, ListCdagC_1, updated,signvalue,activeflavor) 455 456 !Arguments ------------------------------------ 457 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 458 TYPE(MatrixHyb) , INTENT(IN ) :: Mmatrix 459 TYPE(ListCdagC) , INTENT(IN ) :: ListCdagC_1(op%nflavors) 460 DOUBLE PRECISION , INTENT(IN ) :: signvalue 461 LOGICAL , INTENT(IN ) :: updated 462 INTEGER, OPTIONAL , INTENT(IN ) :: activeflavor 463 !Local variables ------------------------------ 464 INTEGER :: iC 465 INTEGER :: iCdag 466 INTEGER :: tail 467 INTEGER :: tailbis 468 !INTEGER :: index 469 INTEGER :: idx_old 470 INTEGER :: old_size 471 ! INTEGER :: omegaSamples 472 ! INTEGER :: iomega 473 INTEGER :: iflavor 474 INTEGER :: iflavorbis 475 INTEGER :: iC_m,iC_m_add 476 INTEGER :: iCdag_m,iCdag_m_add 477 INTEGER :: stail !,ii 478 ! DOUBLE PRECISION :: pi_invBeta 479 DOUBLE PRECISION :: mbeta_two 480 DOUBLE PRECISION :: beta 481 DOUBLE PRECISION :: beta_tc 482 DOUBLE PRECISION :: tcbeta_tc 483 DOUBLE PRECISION :: inv_dt 484 DOUBLE PRECISION :: tC,tc_phys 485 DOUBLE PRECISION :: tCdag 486 DOUBLE PRECISION :: time 487 DOUBLE PRECISION :: signe,signe2 488 DOUBLE PRECISION :: argument 489 INTEGER :: iflavorbegin,iflavorend,prtopt 490 !DOUBLE PRECISION :: taupi_invbeta 491 !COMPLEX(KIND=8) :: cargument 492 !COMPLEX(2*8) :: base_exp 493 !COMPLEX(2*8) :: increm_exp 494 !write(6,*) "measHybrid" 495 prtopt=0 496 IF ( op%set .EQV. .FALSE. ) & 497 CALL ERROR("GreenHyboffdiag_measHybrid : green operator not set ") 498 stail=0 499 do iflavor=1,op%nflavors 500 stail=stail + ListCdagC_1(iflavor)%tail 501 enddo 502 iflavorbegin = 1 503 iflavorend = op%nflavors 504 505 if(present(activeflavor)) then 506 if(activeflavor.ne.0) then 507 !sui!write(6,*) "measHybrid activeflavor",activeflavor 508 iflavorbegin = activeflavor 509 iflavorend = activeflavor 510 endif 511 endif 512 513 IF ( stail .NE. Mmatrix%tail ) & 514 CALL ERROR("GreenHyboffdiag_measHybrid : ListCdagC & M unconsistent ") 515 516 IF ( updated .EQV. .TRUE. ) THEN ! NEW change in the configuration 517 ! FIXME SHOULD be much more faster 518 519 520 ! write(6,*) "LKLLL2b" 521 SELECT CASE(op%iTech) 522 CASE (GREENHYB_TAU) 523 argument = DBLE(op%factor) 524 ! At the beginning old_size=0, then it increases 525 ! until 526 ! for all values of iC, increment green%oper with the value of the 527 ! Green's function in listDBLE(iC) obtained from previous iteration 528 ! (below) 529 ! =============================================================== 530 ! An update has been done. So the Green's function will change 531 ! It is thus the good moment to store the previous Green's 532 ! function with argument, the number of times this Green's 533 ! function has been constant 534 ! =============================================================== 535 DO iflavor=1, op%nflavors 536 DO iflavorbis=1, op%nflavors 537 old_size = op%map(iflavor,iflavorbis)%tail 538 !write(6,*) "size listDBLE",size(op%map(iflavor,iflavorbis)%listDBLE) 539 !sui!write(6,*) " measHybrid",old_size,iflavor,iflavorbis 540 DO iC = 1, old_size 541 if(op%map(iflavor,iflavorbis)%listINT(iC)==0) then 542 !write(6,*) "listINT(iC)=",iC,op%map(iflavor,iflavorbis)%listINT(iC) 543 endif 544 !write(6,*) " measHybrid iflavor,iflavorbis,iC listINT ",iflavor,iflavorbis,iC,op%map(iflavor,iflavorbis)%listINT(iC) 545 !write(6,*) " measHybrid listDBLE ",iflavor,iflavorbis,iC,op%map(iflavor,iflavorbis)%listDBLE(iC),argument 546 op%oper(op%map(iflavor,iflavorbis)%listINT(iC),iflavor,iflavorbis) = & 547 op%oper(op%map(iflavor,iflavorbis)%listINT(iC),iflavor,iflavorbis) & 548 + op%map(iflavor,iflavorbis)%listDBLE(iC) * op%signvalueold * argument 549 !if(op%map(iflavor,iflavorbis)%listINT(iC)==1.and.iflavor==iflavorbis) then 550 ! if(iflavor==iflavorbis) then 551 ! !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 552 ! endif 553 ! if(iflavor==1.and.iflavorbis==6.and.op%map(iflavor,iflavorbis)%listINT(iC)==1) then 554 ! !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 555 ! endif 556 ! if(iflavor==6.and.iflavorbis==1.and.op%map(iflavor,iflavorbis)%listINT(iC)==1) then 557 ! !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 558 ! endif 559 END DO 560 ! tail**2 is the number of possible t-t' 561 ! MapHyb_setSize with resize map tail*tail will thus be the new 562 ! op%map%tail 563 ! update size of map and map%tail 564 CALL MapHyb_setSize(op%map(iflavor,iflavorbis),& 565 & ListCdagC_1(iflavor)%tail*ListCdagC_1(iflavorbis)%tail) 566 END DO 567 END DO 568 op%signvaluemeas = op%signvaluemeas + op%signvalueold * argument 569 op%measurements = op%measurements + op%factor 570 !sui!write(6,*) " measurements", op%measurements 571 !sui! write(6,*) " signvaluemeas",op%signvaluemeas,op%signvalueold*argument 572 !sui! write(6,*) " signvaluemeas/measurements",op%signvaluemeas/op%measurements 573 574 ! This is new measurement, thus op%factor should be put to one 575 op%factor = 1 576 ! write(6,*) "LKLLL2C" 577 578 579 ! initialized index idx_old for the doubles loops over flavors and segments. 580 581 ! setup usefull quantities 582 beta = op%beta 583 mbeta_two = -(beta*0.5d0) 584 inv_dt = op%inv_dt 585 586 ! WARNING time is not the time but just a temporary variable. 587 ! Index Time has been calculated previously and is in mat_tau 588 589 ! initialized index for each annihilation time of a segment for a given flavor 590 591 ! initialized index for each creation time of a segment for another flavor 592 593 iC_m=0 594 iC_m_add=0 595 DO iflavor=1,op%nflavors 596 tail=ListCdagC_1(iflavor)%tail 597 !write(6,*) " measHybrid iflavor",iflavor,tail 598 599 iCdag_m=0 600 iCdag_m_add=0 601 DO iflavorbis=1,op%nflavors 602 tailbis=ListCdagC_1(iflavorbis)%tail 603 !write(6,*) " measHybrid iflavorbis",iflavorbis,tailbis 604 idx_old = 0 605 606 DO iC = 1, tail 607 ! tC is the annihilation (C_) time for segment iC and flavor iflavor 608 !------------------------------------------------------------------- 609 tC = ListCdagC_1(iflavor)%list(iC,C_) 610 611 !iC_m=iC_m+1 ! for Mmatrix%mat 612 ! For each flavor iflavor, iC start at \sum_{iflavor1<iflavor} tail(iflavor1) 613 ! It thus explains the presence of iC_m_add (same below for iCdag_m_add) 614 ! --------------------------------------------------------------------------------- 615 iC_m=iC_m_add+iC 616 beta_tc = beta - tC 617 tcbeta_tc = tC * beta_tc 618 619 !write(6,*) " measHybrid iC_m",iC_m 620 !write(6,*) " measHybrid tailbis",tailbis 621 DO iCdag = 1, tailbis 622 !iCdag_m=iCdag_m+1 623 iCdag_m=iCdag_m_add+iCdag 624 !write(6,*) " measHybrid iCdag_m",iCdag_m 625 626 ! tCdag is the creation time for segment iCdag and flavor iflavorbis 627 tCdag = ListCdagC_1(iflavorbis)%list(iCdag,Cdag_) 628 629 ! --- time is equivalent to time=(tc-tcdag)*(beta-tc) and is only 630 ! --- useful for signe 631 time = tcbeta_tc - tCdag*beta_tc 632 633 !signe = SIGN(1.d0,time) 634 !time = time + (signe-1.d0)*mbeta_two 635 !signe = signe * SIGN(1.d0,beta-tC) 636 !signe = SIGN(1.d0,time) * SIGN(1.d0,beta-tC) 637 tc_phys=tc 638 if(tc>beta) tc_phys=tc-beta 639 signe2=SIGN(1.d0,tc_phys-tcdag) 640 641 if(iflavor==iflavorbis) signe = SIGN(1.d0,time) 642 if(iflavor/=iflavorbis) signe = signe2 643 ! signe = SIGN(1.d0,tc-tcdag) 644 ! --- tc>tcdag and beta>tc signe=1 ! segment in the middle or antisegment at the edge 645 ! ! tc-tcdag > 0 646 ! --- tc<tcdag and beta<tc signe=1 ! never 647 ! --- tc>tcdag and beta<tc signe=-1 ! segment at the edges 648 ! ! tc'-tcdag < 0 (with tc'=tc-beta) -> signe < 0 649 ! --- tc<tcdag and beta>tc signe=-1 ! antisegment in the middle 650 ! ! tc-tcdag < 0 (with tc'=tc-beta) -> signe < 0 651 ! 22/09/14: 652 ! ListCdagC_1 is the list of segment, so we are dealing 653 ! only with segment here. However all combination of Cdag 654 ! and C are taken, this it is possible that tc<tcdag 655 656 ! 21/10/14: Wagt is important are the true times (between 657 ! 0 and beta). If tauC>tauCdag signe=+1 658 ! If tauC<tauCdag signe=-1 659 ! if(tc<tcdag.and.(iflavor==iflavorbis)) then 660 ! write(6,*) ListCdagC_1(iflavorbis)%tail 661 ! do ii=1, ListCdagC_1(iflavorbis)%tail 662 ! write(6,*) ListCdagC_1(iflavorbis)%list(ii,1), ListCdagC_1(iflavorbis)%list(ii,2) 663 ! enddo 664 ! write(6,*) "tc<tcdag", tc,tcdag,beta,iflavor,iflavorbis 665 ! stop 666 ! endif 667 668 if(tc-tcdag>beta) then 669 ! write(6,*) " tc-tcdag > beta ", tcdag-tc,beta 670 endif 671 !if(tc>beta) then 672 ! write(6,*) " TC>BETA" 673 ! write(6,*) " iflavor,iflavorbis",iflavor,iflavorbis 674 ! write(6,*) " ic,icdag ",ic,icdag 675 ! write(6,*) " signe ",signe 676 ! write(6,*) " tc,tcdag ",tc,tcdag 677 ! write(6,*) " Mmatrix%mat ",Mmatrix%mat(iCdag_m,iC_m) 678 ! write(6,*) " Mmatrix%mat_tau ",Mmatrix%mat_tau(iCdag_m,iC_m) 679 !endif 680 ! Si iflavor/=iflavorbis, tc-tcdag can be negative..so in 681 ! this case, on should add beta to tc-tcdag with the minus 682 ! sign. NOT DONE HERE?? 683 684 ! ----- Compute the Green's function as the value of the matrix M for times iCdag and iC. 685 argument = signe*Mmatrix%mat(iCdag_m,iC_m) 686 687 !index = INT( ( time * inv_dt ) + 1.5d0 ) 688 !IF (index .NE. Mmatrix%mat_tau(iCdag,iC)) THEN 689 ! WRITE(*,*) index, Mmatrix%mat_tau(iCdag,iC) 690 !! CALL ERROR("Plantage") 691 !END IF 692 693 idx_old = idx_old + 1 694 695 ! --- define the value of listDBLE as a function of idx_old 696 op%map(iflavor,iflavorbis)%listDBLE(idx_old) = argument 697 !write(6,*) " measHybrid listDBLE2 ",iflavor,iflavorbis,idx_old,argument 698 !op%map%listINT(idx_old) = index 699 700 ! --- define the new corresponding value of listINT(idx_old) from mat_tau (integers) 701 ! --- idx_old has no meaning but listINT(idx_old) has. 702 op%map(iflavor,iflavorbis)%listINT(idx_old) = Mmatrix%mat_tau(iCdag_m,iC_m) 703 !write(6,*) " measHybrid idx_old listINT ",idx_old,op%map(iflavor,iflavorbis)%listINT(idx_old) 704 !write(6,*) " measHybrid iCdag_m, iC_m, mat_tau",iCdag_m,iC_m,Mmatrix%mat_tau(iCdag_m,iC_m) 705 ! if(iflavor==1.and.iflavorbis==2.and.op%map(iflavor,iflavorbis)%listINT(idx_old)==1) then 706 ! !prt!if(prtopt==1) write(6,*) "---------------------------" 707 ! !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 708 ! !prt!if(prtopt==1) write(6,*) " ", tc-tcdag,tc_phys-tcdag 709 ! do ii=1, tail 710 ! !prt!if(prtopt==1) write(6,*) ii, ListCdagC_1(iflavor)%list(ii,1), ListCdagC_1(iflavor)%list(ii,2) 711 ! enddo 712 ! do ii=1, tailbis 713 ! !prt!if(prtopt==1) write(6,*) ii, ListCdagC_1(iflavorbis)%list(ii,1), ListCdagC_1(iflavorbis)%list(ii,2) 714 ! enddo 715 ! !prt!if(prtopt==1) write(6,*) "---------------------------" 716 ! endif 717 ! if(iflavor==1.and.iflavorbis==2.and.op%map(iflavor,iflavorbis)%listINT(idx_old)==999) then 718 ! !prt!if(prtopt==1) write(66,*) "---------------------------" 719 ! !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 720 ! !prt!if(prtopt==1) write(66,*) " ", tc-tcdag,tc_phys-tcdag 721 ! do ii=1, tail 722 ! !prt!if(prtopt==1) write(66,*) ii, ListCdagC_1(iflavor)%list(ii,1), ListCdagC_1(iflavor)%list(ii,2) 723 ! enddo 724 ! do ii=1, tailbis 725 ! !prt!if(prtopt==1) write(66,*) ii, ListCdagC_1(iflavorbis)%list(ii,1), ListCdagC_1(iflavorbis)%list(ii,2) 726 ! enddo 727 ! !prt!if(prtopt==1) write(66,*) "---------------------------" 728 ! endif 729 ! !if(iflavor==2.and.iflavorbis==1.and.op%map(iflavor,iflavorbis)%listINT(idx_old)==1) then 730 ! !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 731 !endif 732 733 734 END DO 735 END DO 736 ! do ii=1,tail*tailbis 737 ! !write(6,*) " measHybrid ii,op%map(iflavor,iflavorbis)%listINT(ii)", ii,op%map(iflavor,iflavorbis)%listINT(ii) 738 ! enddo 739 iCdag_m_add=iCdag_m_add+tailbis 740 END DO ! iflavorbis 741 iC_m_add=iC_m_add+tail 742 END DO ! iflavor 743 op%signvalueold = signvalue 744 ! write(6,*) "LKLLL2D" 745 CASE (GREENHYB_OMEGA) 746 ! argument = DBLE(op%factor) 747 ! DO iomega = 1, omegaSamples 748 ! op%oper_w(iomega) = op%oper_w(iomega) + op%oper_w_old(iomega) * argument 749 ! END DO 750 ! op%measurements = op%measurements + op%factor 751 752 ! op%factor = 1 753 ! beta = op%beta 754 ! mbeta_two = -(beta*0.5d0) 755 ! pi_invBeta = ACOS(-1.d0)/beta 756 ! omegaSamples = op%samples-1 757 ! DO iC = 1, tail 758 ! tC = ListCdagC_1%list(iC,C_) 759 ! DO iCdag = 1, tail 760 ! tCdag = ListCdagC_1%list(iCdag,Cdag_) 761 ! time = tC - tCdag 762 763 ! signe = SIGN(1.d0,time) 764 ! time = time + (signe-1.d0)*mbeta_two 765 ! signe = signe * SIGN(1.d0,beta-tC) 766 ! argument = signe*Mmatrix%mat(iCdag,iC) 767 768 ! DO iomega = 1, omegaSamples 769 ! !op%oper_w_old(iomega) = Mmatrix%mat_tau(iCdag,iC)*CMPLX(0.d0,argument) 770 ! op%oper_w_old(iomega) = EXP(CMPLX(0.d0,op%omega(iomega)*time))*CMPLX(0.d0,argument) 771 ! END DO 772 ! END DO 773 ! END DO 774 END SELECT 775 ELSE 776 op%factor = op%factor + 1 777 END IF 778 END SUBROUTINE GreenHyboffdiag_measHybrid
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_print [ Functions ]
NAME
GreenHyboffdiag_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
op=Green ostream=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1734 SUBROUTINE GreenHyboffdiag_print(op, ostream) 1735 1736 !Arguments ------------------------------------ 1737 TYPE(GreenHyboffdiag), INTENT(IN) :: op 1738 INTEGER, OPTIONAL , INTENT(IN) :: ostream 1739 !Local variables ------------------------------ 1740 INTEGER :: ostream_val 1741 INTEGER :: isample 1742 INTEGER :: samples 1743 INTEGER :: iflavor1 1744 INTEGER :: iflavor2 1745 1746 1747 IF ( op%set .EQV. .FALSE. ) & 1748 CALL ERROR("GreenHyboffdiag_print : green op%operator not set ") 1749 1750 IF ( PRESENT(ostream) ) THEN 1751 ostream_val = ostream 1752 ELSE 1753 ostream_val = 66 1754 OPEN(UNIT=ostream_val,FILE="Green.dat") 1755 END IF 1756 1757 samples = op%samples 1758 1759 DO iflavor1=1,op%nflavors 1760 DO iflavor2=1,op%nflavors 1761 WRITE(ostream_val,'(a,i3,a,i3,a)') "## (iflavor1,iflavor2)= (", iflavor1,",",iflavor2,")" 1762 DO isample = 1, samples 1763 WRITE(ostream_val,*) DBLE(isample-1)*op%delta_t, op%oper(isample,iflavor1,iflavor2) 1764 END DO 1765 WRITE(ostream_val,*) 1766 END DO 1767 END DO 1768 1769 IF ( .NOT. PRESENT(ostream) ) & 1770 CLOSE(ostream_val) 1771 END SUBROUTINE GreenHyboffdiag_print
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_reset [ Functions ]
NAME
GreenHyboffdiag_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
op=Green
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
311 SUBROUTINE GreenHyboffdiag_reset(op) 312 313 !Arguments ------------------------------------ 314 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 315 316 CALL GreenHyboffdiag_clear(op) 317 op%setMk = 0 318 op%Mk = 0.d0 319 op%setT = .FALSE. 320 op%setW = .FALSE. 321 END SUBROUTINE GreenHyboffdiag_reset
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setMoments [ Functions ]
NAME
GreenHyboffdiag_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
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
CHI Will be filled automatically by the parent script
SOURCE
987 SUBROUTINE GreenHyboffdiag_setMoments(op,iflavor1,iflavor1b,u1,u2,u3) 988 989 !Arguments ------------------------------------ 990 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 991 DOUBLE PRECISION, INTENT(IN ) :: u1 992 DOUBLE PRECISION, INTENT(IN ) :: u2 993 DOUBLE PRECISION, INTENT(IN ) :: u3 994 INTEGER , INTENT(IN ) :: iflavor1 995 INTEGER , INTENT(IN ) :: iflavor1b 996 997 if(iflavor1==iflavor1b) then 998 op%Mk(iflavor1,iflavor1b,1) = -1.d0 999 ! c_a(3)=-d1-mu*mu-2(-mu)(\sum_{b.ne.a} Uab nb) 1000 op%Mk(iflavor1,iflavor1b,3) = op%Mk(iflavor1,iflavor1b,3) - 2.d0*(op%Mk(iflavor1,iflavor1b,2)*u1) 1001 1002 ! c_a(2)=-mu+\sum_{b.ne.a} Uab n_b 1003 op%Mk(iflavor1,iflavor1b,2) = op%Mk(iflavor1,iflavor1b,2) + u1 1004 !sui!write(6,*) "setmiments",iflavor1,iflavor1b,u1 1005 1006 ! c_a(3)=c_a(3) + \sum Uab^2 nb + \sum Uba Uca <nbnc> 1007 ! ie c_a(3)=-d1+mu*mu-2mu*\sumb Uab nb + \sum Uab^2 nb + \sum Uba Uca <nbnc> 1008 op%Mk(iflavor1,iflavor1b,3) = op%Mk(iflavor1,iflavor1b,3) - u2 1009 else 1010 op%Mk(iflavor1,iflavor1b,1) = 0.d0 1011 op%Mk(iflavor1,iflavor1b,2) = u3 1012 op%Mk(iflavor1,iflavor1b,3) = 0.d0 1013 endif 1014 !write(6,*) "mom",iflavor1,iflavor1b, op%Mk(iflavor1,iflavor1b,:) 1015 1016 op%setMk = op%setMk + 1 1017 1018 END SUBROUTINE GreenHyboffdiag_setMoments
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setMuD1 [ Functions ]
NAME
GreenHyboffdiag_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
op=Green mu=energy level (irrespectige with fermi level) d1=first moment of hybridization function ("K")
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
923 SUBROUTINE GreenHyboffdiag_setMuD1(op,iflavor,iflavor2,mu,d1) 924 925 !Arguments ------------------------------------ 926 !Arguments ------------------------------------ 927 !scalars 928 DOUBLE PRECISION, INTENT(IN ) :: mu 929 DOUBLE PRECISION, INTENT(IN ) :: d1 930 INTEGER , INTENT(IN ) :: iflavor 931 INTEGER , INTENT(IN ) :: iflavor2 932 !type 933 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 934 !Local variables ------------------------------------ 935 DOUBLE PRECISION :: mu2 936 !********************************************************************* 937 938 ABI_UNUSED((/d1/)) 939 940 mu2=0 941 if(iflavor==iflavor2) mu2=mu 942 943 if(iflavor==iflavor2) then 944 op%Mk(iflavor,iflavor2,3) = -d1-(mu*mu) 945 !op%Mk(iflavor,iflavor2,3) = -(mu*mu) 946 op%Mk(iflavor,iflavor2,2) = -mu 947 !sui!write(6,*) "setmud1",iflavor,iflavor2, op%Mk(iflavor,iflavor2,2), op%Mk(iflavor,iflavor2,3) 948 else 949 op%Mk(iflavor,iflavor2,3) = 0.d0 950 op%Mk(iflavor,iflavor2,2) = 0.d0 951 endif 952 op%setMk = op%setMk + 1 953 !write(6,*) "mom1",op%Mk(iflavor,iflavor2,:) 954 END SUBROUTINE GreenHyboffdiag_setMuD1
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setN [ Functions ]
NAME
GreenHyboffdiag_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
op=Green N=number of electrons
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
857 SUBROUTINE GreenHyboffdiag_setN(op,N) 858 859 !Arguments ------------------------------------ 860 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 861 DOUBLE PRECISION , INTENT(IN ) :: N(op%nflavors) 862 INTEGER :: iflavor,iflavor2 863 !COMPLEX(KIND=8) :: tmpoper 864 865 IF ( op%set .EQV. .FALSE. ) & 866 CALL ERROR("GreenHyboffdiag_setN: green op%operator not set ") 867 DO iflavor=1, op%nflavors 868 ! write(6,*) "iflavor",-N(iflavor)*op%signvaluemeas ,2*op%oper(op%samples,iflavor,iflavor),(N(iflavor)-1.d0)*op%signvaluemeas 869 ! the mulplication by signvaluemeas is necessary because N is 870 ! exactly the number of electrons in the flavor iflavor whereas 871 ! op%oper is not exact, because it still has to be divided by 872 ! signvaluemeas after the MPIREDUCE 873 op%oper(1,iflavor,iflavor) = (N(iflavor) - 1.d0)*op%signvaluemeas 874 op%oper(op%samples,iflavor,iflavor) = - N(iflavor)*op%signvaluemeas 875 !op%oper(op%samples,iflavor,iflavor) = 2*op%oper(op%samples,iflavor,iflavor) 876 !op%oper(1,iflavor,iflavor) = 2*op%oper(1,iflavor,iflavor) 877 DO iflavor2=1, op%nflavors 878 if(iflavor/=iflavor2) then 879 ! UNEXPLAINED but MANDATORY to have exact results for U=0 nspinor=4 with pawspnorb=0 880 ! Correction: The fact 2 is necessary for edge points because the points are at the 881 ! edges. 882 ! It is of course necessary to fulfill exact results (U=0). 883 !tmpoper=(op%oper(op%samples,iflavor,iflavor2)-op%oper(1,iflavor,iflavor2)) 884 !op%oper(op%samples,iflavor,iflavor2) = tmpoper 885 !op%oper(1,iflavor,iflavor2) = -tmpoper 886 !op%oper(op%samples,iflavor,iflavor2) = 2*op%oper(op%samples,iflavor,iflavor2) 887 !op%oper(1,iflavor,iflavor2) = 2*op%oper(1,iflavor,iflavor2) 888 op%oper(op%samples,iflavor,iflavor2) = 2*op%oper(op%samples,iflavor,iflavor2) 889 op%oper(1,iflavor,iflavor2) = 2*op%oper(1,iflavor,iflavor2) 890 endif 891 ENDDO 892 ENDDO 893 END SUBROUTINE GreenHyboffdiag_setN
ABINIT/m_GreenHyboffdiag/GreenHyboffdiag_setOperW [ Functions ]
NAME
GreenHyboffdiag_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
op=Green Gomega=Input values
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
402 SUBROUTINE GreenHyboffdiag_setOperW(op, Gomega) 403 404 !Arguments ------------------------------------ 405 TYPE(GreenHyboffdiag) , INTENT(INOUT) :: op 406 COMPLEX(KIND=8), DIMENSION(:,:,:), INTENT(IN ) :: Gomega 407 !Loval variables ------------------------------ 408 INTEGER :: tail 409 410 tail = SIZE(Gomega,1) 411 IF ( .NOT. op%set ) & 412 CALL ERROR("GreenHyboffdiag_setOperW : Uninitialized GreenHyboffdiag structure") 413 IF ( ALLOCATED(op%oper_w) ) THEN 414 IF ( SIZE(op%oper_w) .NE. tail ) THEN 415 FREE(op%oper_w) 416 MALLOC(op%oper_w,(1:tail,op%nflavors,op%nflavors)) 417 END IF 418 ELSE 419 MALLOC(op%oper_w,(1:tail,op%nflavors,op%nflavors)) 420 END IF 421 op%oper_w(:,:,:) = Gomega(:,:,:) 422 op%Wmax = tail 423 op%setW = .TRUE. 424 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-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
69 TYPE GreenHyboffdiag 70 71 LOGICAL :: set = .FALSE. 72 ! True if variable of type GreenHyboffdiag is initialized 73 74 LOGICAL :: setT = .FALSE. 75 ! True if variable oper contains data 76 77 LOGICAL :: setW = .FALSE. 78 ! True if variable oper_w contains data 79 80 LOGICAL :: have_MPI = .FALSE. 81 ! True if MPI is used. 82 83 INTEGER :: setMk = 0 84 ! setMk=0 is moments for Fourier transform are not computed 85 86 INTEGER :: samples 87 ! samples=imaginary time slices (dmftqmc_l+1) 88 89 INTEGER :: measurements 90 ! number of measurements for the Green's function 91 92 INTEGER :: factor 93 ! if the move is not accepted, the statistic weight has to be 94 ! increased for the current configuration. 95 96 INTEGER :: MY_COMM 97 ! MPI Communicator 98 99 INTEGER :: size 100 ! size=1 101 102 INTEGER :: rank 103 ! rank=0 104 105 INTEGER :: Wmax 106 ! samples-1 if frequency Green's function 107 108 INTEGER :: iTech 109 ! Precise if Frequency Green's function is computed or not 110 111 INTEGER :: nflavors 112 ! Number of flavors 113 114 DOUBLE PRECISION :: beta 115 ! Inverse of temperature 116 117 DOUBLE PRECISION :: inv_beta 118 ! Temperature 119 120 DOUBLE PRECISION :: delta_t 121 ! 1/inv_dt 122 123 DOUBLE PRECISION :: inv_dt 124 ! (samples-1)/beta 125 DOUBLE PRECISION :: signvaluemeas 126 127 DOUBLE PRECISION :: signvalueold 128 129 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: oper 130 ! oper(samples) 131 132 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: omega 133 ! omega(Wmax) 134 135 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: Mk 136 ! Moments for FT 137 138 COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:,:,:) :: oper_w 139 ! Frequency Green's function 140 141 COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:) :: oper_w_old 142 ! Old frequency Green's function (not used) 143 144 TYPE(Vector) :: oper_old 145 ! useless data 146 147 TYPE(VectorInt) :: index_old 148 ! useless data 149 150 TYPE(MapHyb), ALLOCATABLE, DIMENSION(:,:) :: map 151 ! value of time and Green's functions computed in GreenHyboffdiag_measHybrid 152 ! These values are used to fill op%oper in the same routine. 153 154 END TYPE GreenHyboffdiag