TABLE OF CONTENTS
- ABINIT/m_Ctqmc
- ABINIT/m_Ctqmc/Ctqmc_allocateAll
- ABINIT/m_Ctqmc/Ctqmc_allocateOpt
- ABINIT/m_Ctqmc/Ctqmc_clear
- ABINIT/m_Ctqmc/Ctqmc_computeF
- ABINIT/m_Ctqmc/Ctqmc_destroy
- ABINIT/m_Ctqmc/Ctqmc_getD
- ABINIT/m_Ctqmc/Ctqmc_getE
- ABINIT/m_Ctqmc/Ctqmc_getGreen
- ABINIT/m_Ctqmc/Ctqmc_getResult
- ABINIT/m_Ctqmc/Ctqmc_init
- ABINIT/m_Ctqmc/Ctqmc_loop
- ABINIT/m_Ctqmc/Ctqmc_measCorrelation
- ABINIT/m_Ctqmc/Ctqmc_measN
- ABINIT/m_Ctqmc/Ctqmc_measPerturbation
- ABINIT/m_Ctqmc/Ctqmc_printAll
- ABINIT/m_Ctqmc/Ctqmc_printCorrelation
- ABINIT/m_Ctqmc/Ctqmc_printD
- ABINIT/m_Ctqmc/Ctqmc_printE
- ABINIT/m_Ctqmc/Ctqmc_printGreen
- ABINIT/m_Ctqmc/Ctqmc_printPerturbation
- ABINIT/m_Ctqmc/Ctqmc_printQMC
- ABINIT/m_Ctqmc/Ctqmc_printSpectra
- ABINIT/m_Ctqmc/Ctqmc_reset
- ABINIT/m_Ctqmc/Ctqmc_run
- ABINIT/m_Ctqmc/Ctqmc_setG0wTab
- ABINIT/m_Ctqmc/Ctqmc_setMu
- ABINIT/m_Ctqmc/Ctqmc_setParameters
- ABINIT/m_Ctqmc/Ctqmc_setSeed
- ABINIT/m_Ctqmc/Ctqmc_setSweeps
- ABINIT/m_Ctqmc/Ctqmc_setU
- ABINIT/m_Ctqmc/Ctqmc_symmetrizeGreen
- ABINIT/m_Ctqmc/Ctqmc_tryAddRemove
- ABINIT/m_Ctqmc/Ctqmc_trySwap
- m_Ctqmc/Ctqmc
ABINIT/m_Ctqmc [ Modules ]
NAME
m_Ctqmc
FUNCTION
Manage and drive all the CTQMC Should not be used if you don't know what you do Please use CtqmcInterface
COPYRIGHT
Copyright (C) 2013-2022 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
25 #include "defs.h" 26 27 MODULE m_Ctqmc 28 29 USE m_Global 30 USE m_GreenHyb 31 USE m_BathOperator 32 USE m_ImpurityOperator 33 USE m_Stat 34 USE m_FFTHyb 35 USE m_OurRng 36 USE m_Vector 37 #ifdef HAVE_MPI2 38 USE mpi 39 #endif 40 41 IMPLICIT NONE
ABINIT/m_Ctqmc/Ctqmc_allocateAll [ Functions ]
NAME
Ctqmc_allocateAll
FUNCTION
Allocate all non option varibales
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
622 SUBROUTINE Ctqmc_allocateAll(this) 623 624 !Arguments ------------------------------------ 625 TYPE(Ctqmc), INTENT(INOUT) :: this 626 !Local variables ------------------------------ 627 INTEGER :: flavors 628 629 IF ( .NOT. this%para ) & 630 CALL ERROR("Ctqmc_allocateAll : Ctqmc_setParameters never called ") 631 632 flavors = this%flavors 633 634 DT_FREEIF(this%Greens) 635 DT_MALLOC(this%Greens,(1:flavors)) 636 637 FREEIF(this%measN) 638 MALLOC(this%measN,(1:4,1:flavors)) 639 this%measN = 0.d0 640 641 FREEIF(this%measDE) 642 MALLOC(this%measDE,(1:flavors,1:flavors) ) 643 this%measDE = 0.d0 644 645 FREEIF(this%mu) 646 MALLOC(this%mu,(1:flavors) ) 647 this%mu = 0.d0 648 END SUBROUTINE Ctqmc_allocateAll
ABINIT/m_Ctqmc/Ctqmc_allocateOpt [ Functions ]
NAME
Ctqmc_allocateOpt
FUNCTION
allocate all option variables
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
676 SUBROUTINE Ctqmc_allocateOpt(this) 677 678 !Arguments ------------------------------------ 679 TYPE(Ctqmc), INTENT(INOUT) :: this 680 !Local variables ------------------------------ 681 INTEGER :: i 682 INTEGER :: j 683 INTEGER :: k 684 685 IF ( .NOT. this%para ) & 686 CALL ERROR("Ctqmc_allocateOpt : Ctqmc_setParameters never called ") 687 688 IF ( this%opt_analysis .EQ. 1 ) THEN 689 FREEIF(this%measCorrelation) 690 MALLOC(this%measCorrelation,(1:this%samples+1,1:3,1:this%flavors)) 691 this%measCorrelation = 0.d0 692 END IF 693 694 IF ( this%opt_order .GT. 0 ) THEN 695 FREEIF(this%measPerturbation) 696 MALLOC(this%measPerturbation,(1:this%opt_order,1:this%flavors)) 697 this%measPerturbation = 0.d0 698 END IF 699 700 IF ( this%opt_histo .GT. 0 ) THEN 701 FREEIF(this%occup_histo_time) 702 MALLOC(this%occup_histo_time,(1:this%flavors+1)) 703 this%occup_histo_time= 0.d0 704 END IF 705 706 IF ( this%opt_noise .EQ. 1 ) THEN 707 IF ( ALLOCATED(this%measNoiseG) ) THEN 708 DO i=1,2 709 DO j = 1, this%flavors 710 DO k= 1, this%samples+1 711 CALL Vector_destroy(this%measNoiseG(k,j,i)) 712 END DO 713 END DO 714 END DO 715 DT_FREE(this%measNoiseG) 716 END IF 717 DT_MALLOC(this%measNoiseG,(1:this%samples+1,1:this%flavors,1:2)) 718 !DO i=1,2 719 DO j = 1, this%flavors 720 DO k= 1, this%samples+1 721 CALL Vector_init(this%measNoiseG(k,j,1),CTQMC_SLICE1) 722 END DO 723 END DO 724 DO j = 1, this%flavors 725 DO k= 1, this%samples+1 726 CALL Vector_init(this%measNoiseG(k,j,2),CTQMC_SLICE1*CTQMC_SLICE2+1) ! +1 pour etre remplacer ceil 727 END DO 728 END DO 729 !END DO 730 FREEIF(this%abNoiseG) 731 MALLOC(this%aBNoiseG,(1:2,1:this%samples+1,this%flavors)) 732 this%abNoiseG = 0.d0 733 END IF 734 735 IF (this%opt_spectra .GE. 1 ) THEN 736 FREEIF(this%density) 737 !MALLOC(this%density,(1:this%thermalization,1:this%flavors)) 738 i = CEILING(DBLE(this%thermalization+this%sweeps)/DBLE(this%measurements*this%opt_spectra)) 739 MALLOC(this%density,(1:this%flavors+1,1:i)) 740 this%density = 0.d0 741 END IF 742 !#endif 743 END SUBROUTINE Ctqmc_allocateOpt
ABINIT/m_Ctqmc/Ctqmc_clear [ Functions ]
NAME
Ctqmc_clear
FUNCTION
clear a ctqmc run
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
937 SUBROUTINE Ctqmc_clear(this) 938 939 !Arguments ------------------------------------ 940 TYPE(Ctqmc), INTENT(INOUT) :: this 941 !Local variables ------------------------------ 942 INTEGER :: i 943 INTEGER :: j 944 INTEGER :: k 945 946 this%measN(1,:) = 0.d0 947 this%measN(2,:) = 0.d0 948 !Do not set measN(3,:) to 0 to avoid erasing N between therm and ctqmc 949 this%measN(4,:) = 0.d0 950 this%measDE = 0.d0 951 ! this%seg_added = 0.d0 952 ! this%anti_added = 0.d0 953 ! this%seg_removed = 0.d0 954 ! this%anti_removed = 0.d0 955 ! this%seg_sign = 0.d0 956 ! this%anti_sign = 0.d0 957 this%stats(:) = 0.d0 958 this%swap = 0.d0 959 this%runTime = 0.d0 960 this%modGlobalMove(2) = 0 961 CALL Vector_clear(this%measNoise(1)) 962 CALL Vector_clear(this%measNoise(2)) 963 !#ifdef CTCtqmc_CHECK 964 this%errorImpurity = 0.d0 965 this%errorBath = 0.d0 966 !#endif 967 DO j = 1, this%flavors 968 CALL GreenHyb_clear(this%Greens(j)) 969 END DO 970 !#ifdef CTCtqmc_ANALYSIS 971 IF ( this%opt_analysis .EQ. 1 .AND. ALLOCATED(this%measCorrelation) ) & 972 this%measCorrelation = 0.d0 973 IF ( this%opt_order .GT. 0 .AND. ALLOCATED(this%measPerturbation) ) & 974 this%measPerturbation = 0.d0 975 IF ( this%opt_noise .EQ. 1 .AND. ALLOCATED(this%measNoiseG) ) THEN 976 DO i=1,2 977 DO j = 1, this%flavors 978 DO k= 1, this%samples+1 979 CALL Vector_clear(this%measNoiseG(k,j,i)) 980 END DO 981 END DO 982 END DO 983 END IF 984 !#endif 985 END SUBROUTINE Ctqmc_clear
ABINIT/m_Ctqmc/Ctqmc_computeF [ Functions ]
NAME
Ctqmc_computeF
FUNCTION
Compute the hybridization function
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc Gomega=G0 to compute F opt_fk=What is Gomega
OUTPUT
F=hybridization function
SIDE EFFECTS
NOTES
SOURCE
1118 SUBROUTINE Ctqmc_computeF(this, Gomega, F, opt_fk) 1119 1120 !Arguments ------------------------------------ 1121 TYPE(Ctqmc) , INTENT(INOUT) :: this 1122 COMPLEX(KIND=8), DIMENSION(:,:), INTENT(IN ) :: Gomega 1123 !INTEGER , INTENT(IN ) :: Wmax 1124 DOUBLE PRECISION, DIMENSION(:,:), INTENT(INOUT) :: F 1125 INTEGER , INTENT(IN ) :: opt_fk 1126 !Local variables ------------------------------ 1127 INTEGER :: flavors 1128 INTEGER :: samples 1129 INTEGER :: iflavor 1130 INTEGER :: iomega 1131 INTEGER :: itau 1132 DOUBLE PRECISION :: pi_invBeta 1133 DOUBLE PRECISION :: K 1134 DOUBLE PRECISION :: re 1135 DOUBLE PRECISION :: im 1136 COMPLEX(KIND=8), DIMENSION(:,:), ALLOCATABLE :: F_omega 1137 TYPE(GreenHyb) :: F_tmp 1138 1139 flavors = this%flavors 1140 1141 samples = this%samples 1142 pi_invBeta = ACOS(-1.d0) / this%beta 1143 this%Wmax=SIZE(Gomega,1) 1144 1145 IF ( this%have_MPI .EQV. .TRUE. ) THEN 1146 CALL GreenHyb_init(F_tmp,samples,this%beta,MY_COMM=this%MY_COMM) 1147 ELSE 1148 CALL GreenHyb_init(F_tmp,samples,this%beta) 1149 END IF 1150 ! K = this%mu 1151 1152 MALLOC(F_omega,(1:this%Wmax,1:flavors)) 1153 1154 !IF ( this%rank .EQ. 0 ) & 1155 !OPEN(UNIT=9876,FILE="K.dat",POSITION="APPEND") 1156 IF ( opt_fk .EQ. 0 ) THEN 1157 DO iflavor = 1, flavors 1158 DO iomega=1,this%Wmax 1159 re = REAL(Gomega(iomega,iflavor)) 1160 im = AIMAG(Gomega(iomega,iflavor)) 1161 F_omega(iomega,iflavor) = CMPLX(-re/(re*re+im*im),im/(re*re+im*im),8) 1162 END DO 1163 END DO 1164 !F_omega = CMPLX(-1.d0,0,8)/Gomega 1165 ELSE 1166 F_omega = Gomega 1167 END IF 1168 1169 DO iflavor = 1, flavors 1170 IF ( this%opt_levels .EQ. 1 ) THEN 1171 K = this%mu(iflavor) 1172 ELSE 1173 K = -REAL(F_omega(this%Wmax, iflavor)) 1174 ! this%mu = K 1175 this%mu(iflavor) = K 1176 END IF 1177 !IF ( this%rank .EQ. 0 ) & 1178 !WRITE(9876,'(I4,2E22.14)') iflavor, K, REAL(-F_omega(this%Wmax, iflavor)) 1179 !IF(this%rank .EQ.0) & 1180 !WRITE(this%ostream,*) "CTQMC K, this%mu = ",K,this%mu(iflavor) 1181 !WRITE(this%ostream,*) "CTQMC beta = ",this%beta 1182 IF ( opt_fk .EQ. 0 ) THEN 1183 DO iomega = 1, this%Wmax 1184 re = REAL(F_omega(iomega,iflavor)) 1185 im = AIMAG(F_omega(iomega,iflavor)) 1186 F_omega(iomega,iflavor) = CMPLX(re + K, im + (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, 8) 1187 !if(iflavor==1.and.this%rank==0) then 1188 !write(224,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(F_omega(iomega,iflavor)),imag(F_omega(iomega,iflavor)) 1189 !write(225,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(Gomega(iomega, iflavor)),imag(Gomega(iomega, iflavor)) 1190 !end if 1191 END DO 1192 ELSE 1193 DO iomega = 1, this%Wmax 1194 F_omega(iomega,iflavor) = F_omega(iomega,iflavor) & 1195 + CMPLX(K, 0.d0, 8) 1196 !if(iflavor==1.and.this%rank==0) then 1197 !write(224,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(F_omega(iomega,iflavor)),imag(F_omega(iomega,iflavor)) 1198 !write(225,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(Gomega(iomega, iflavor)),imag(Gomega(iomega, iflavor)) 1199 !end if 1200 END DO 1201 END IF 1202 K = REAL(CMPLX(0,(2.d0*DBLE(this%Wmax)-1.d0)*pi_invBeta,8)*F_omega(this%Wmax,iflavor)) 1203 CALL GreenHyb_setMuD1(this%Greens(iflavor),this%mu(iflavor),K) 1204 CALL GreenHyb_setOperW(F_tmp,F_omega(:,iflavor)) 1205 !CALL GreenHyb_backFourier(F_tmp,F_omega(:,iflavor)) 1206 CALL GreenHyb_backFourier(F_tmp) 1207 F(1:samples+1,iflavor) = (/ (-F_tmp%oper(samples+1-itau),itau=0,samples) /) 1208 END DO 1209 IF ( this%rank .EQ. 0 ) THEN 1210 open (unit=346,file='Hybridization.dat',status='unknown',form='formatted') 1211 DO iflavor = 1, flavors 1212 write(346,*) "#",iflavor 1213 do itau=1,this%samples+1 1214 write(346,*) itau,F(itau,iflavor) 1215 enddo 1216 write(346,*) 1217 END DO 1218 close(346) 1219 ENDIF 1220 FREE(F_omega) 1221 CALL GreenHyb_destroy(F_tmp) 1222 END SUBROUTINE Ctqmc_computeF
ABINIT/m_Ctqmc/Ctqmc_destroy [ Functions ]
NAME
Ctqmc_destroy
FUNCTION
destroy and deallocate all variables
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3349 SUBROUTINE Ctqmc_destroy(this) 3350 3351 !Arguments ------------------------------------ 3352 TYPE(Ctqmc), INTENT(INOUT) :: this 3353 !Local variables ------------------------------ 3354 INTEGER :: iflavor 3355 INTEGER :: flavors 3356 INTEGER :: i 3357 INTEGER :: j 3358 INTEGER :: k 3359 3360 if ( this%init .EQV. .FALSE. ) RETURN 3361 3362 flavors = this%flavors 3363 3364 CALL ImpurityOperator_destroy(this%Impurity) 3365 CALL BathOperator_destroy(this%Bath) 3366 CALL Vector_destroy(this%measNoise(1)) 3367 CALL Vector_destroy(this%measNoise(2)) 3368 3369 IF ( ALLOCATED(this%Greens) ) THEN 3370 DO iflavor = 1, flavors 3371 CALL GreenHyb_destroy(this%Greens(iflavor)) 3372 END DO 3373 DT_FREE( this%Greens ) 3374 END IF 3375 !#ifdef CTCtqmc_ANALYSIS 3376 FREEIF(this%measCorrelation) 3377 FREEIF(this%measPerturbation) 3378 IF ( this%opt_histo .GT. 0 ) THEN 3379 FREEIF(this%occup_histo_time) 3380 ENDIF 3381 FREEIF(this%measN) 3382 FREEIF(this%measDE) 3383 FREEIF(this%mu) 3384 FREEIF(this%abNoiseG) 3385 IF ( ALLOCATED(this%measNoiseG) ) THEN 3386 DO i=1,2 3387 DO j = 1, this%flavors 3388 DO k= 1, this%samples+1 3389 CALL Vector_destroy(this%measNoiseG(k,j,i)) 3390 END DO 3391 END DO 3392 END DO 3393 DT_FREE(this%measNoiseG) 3394 END IF 3395 FREEIF(this%density) 3396 !#endif 3397 this%ostream = 0 3398 this%istream = 0 3399 3400 this%sweeps = 0 3401 this%thermalization = 0 3402 this%flavors = 0 3403 this%samples = 0 3404 this%beta = 0.d0 3405 ! this%seg_added = 0.d0 3406 ! this%anti_added = 0.d0 3407 ! this%seg_removed = 0.d0 3408 ! this%anti_removed = 0.d0 3409 ! this%seg_sign = 0.d0 3410 ! this%anti_sign = 0.d0 3411 this%stats = 0.d0 3412 this%swap = 0.d0 3413 3414 3415 this%set = .FALSE. 3416 this%done = .FALSE. 3417 this%init = .FALSE. 3418 END SUBROUTINE Ctqmc_destroy
ABINIT/m_Ctqmc/Ctqmc_getD [ Functions ]
NAME
Ctqmc_getD
FUNCTION
get double occupation
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
D=full double occupation
SIDE EFFECTS
NOTES
SOURCE
2710 SUBROUTINE Ctqmc_getD(this, D) 2711 2712 !Arguments ------------------------------------ 2713 TYPE(Ctqmc) , INTENT(IN ) :: this 2714 DOUBLE PRECISION, DIMENSION(:,:), INTENT(OUT) :: D 2715 !Local variables ------------------------------ 2716 INTEGER :: iflavor1 2717 INTEGER :: iflavor2 2718 2719 IF ( SIZE(D,1) .LT. this%flavors .OR. SIZE(D,2) .LT. this%flavors ) & 2720 CALL ERROR("Ctqmc_getD : Dimensions of array D are too small") 2721 2722 D = 0.d0 2723 2724 DO iflavor1 = 1, this%flavors 2725 DO iflavor2 = 1, this%flavors 2726 D(iflavor2,iflavor1) = this%measDE(iflavor2,iflavor1) 2727 END DO 2728 D(iflavor1,iflavor1) = 0.d0 2729 END DO 2730 2731 END SUBROUTINE Ctqmc_getD
ABINIT/m_Ctqmc/Ctqmc_getE [ Functions ]
NAME
Ctqmc_getE
FUNCTION
get interaction energy and noise on it
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
E=interaction energy noise=noise on this value
SIDE EFFECTS
NOTES
SOURCE
2760 SUBROUTINE Ctqmc_getE(this,E,noise) 2761 2762 !Arguments ------------------------------------ 2763 TYPE(Ctqmc) , INTENT(IN ) :: this 2764 DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: E 2765 DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: Noise 2766 2767 IF ( PRESENT(E) ) & 2768 E = this%measDE(1,1) 2769 IF ( PRESENT(Noise) ) & 2770 Noise = SUM(this%Impurity%mat_U)/(this%flavors*(this%flavors-1)) & 2771 * this%a_Noise*(DBLE(this%sweeps)*DBLE(this%size))**this%b_Noise 2772 END SUBROUTINE Ctqmc_getE
ABINIT/m_Ctqmc/Ctqmc_getGreen [ Functions ]
NAME
Ctqmc_getGreen
FUNCTION
Get the full green functions in time and/or frequency
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
Gtau=green function in time Gw=green function in frequency
SIDE EFFECTS
NOTES
SOURCE
2629 SUBROUTINE Ctqmc_getGreen(this, Gtau, Gw) 2630 2631 !Arguments ------------------------------------ 2632 TYPE(Ctqmc) , INTENT(INOUT) :: this 2633 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT(INOUT) :: Gtau 2634 COMPLEX(KIND=8), DIMENSION(:,:), OPTIONAL, INTENT(INOUT) :: Gw 2635 !Local variables ------------------------------ 2636 !INTEGER :: itime 2637 INTEGER :: iflavor1 2638 INTEGER :: iflavor2 2639 INTEGER :: iflavor3 2640 INTEGER :: flavors 2641 DOUBLE PRECISION :: u1 2642 DOUBLE PRECISION :: u2 2643 DOUBLE PRECISION :: Un 2644 DOUBLE PRECISION :: UUnn 2645 2646 flavors = this%flavors 2647 DO iflavor1 = 1, flavors 2648 u1 = 0.d0 2649 u2 = 0.d0 2650 DO iflavor2 = 1, flavors 2651 IF ( iflavor2 .EQ. iflavor1 ) CYCLE 2652 Un = this%Impurity%mat_U(iflavor2,iflavor1) * this%measN(1,iflavor2) 2653 u1 = u1 + Un 2654 u2 = u2 + Un*this%Impurity%mat_U(iflavor2,iflavor1) 2655 DO iflavor3 = 1, flavors 2656 IF ( iflavor3 .EQ. iflavor2 .OR. iflavor3 .EQ. iflavor1 ) CYCLE 2657 UUnn = (this%Impurity%mat_U(iflavor2,iflavor1)*this%Impurity%mat_U(iflavor3,iflavor1)) * this%measDE(iflavor2,iflavor3) 2658 u2 = u2 + UUnn 2659 END DO 2660 END DO 2661 2662 CALL GreenHyb_setMoments(this%Greens(iflavor1),u1,u2) 2663 IF ( PRESENT( Gtau ) ) THEN 2664 Gtau(1:this%samples,iflavor1) = this%Greens(iflavor1)%oper(1:this%samples) 2665 END IF 2666 !write(6,*) "present gw", present(gw) 2667 IF ( PRESENT( Gw ) ) THEN 2668 !write(6,*) "size gw",SIZE(Gw,DIM=2) ,flavors+1 2669 IF ( SIZE(Gw,DIM=2) .EQ. flavors+1 ) THEN 2670 CALL GreenHyb_forFourier(this%Greens(iflavor1), Gomega=Gw(:,iflavor1), omega=Gw(:,this%flavors+1)) 2671 !IF ( this%rank .EQ. 0 ) write(20,*) Gw(:,iflavor1) 2672 ELSE IF ( SIZE(Gw,DIM=2) .EQ. flavors ) THEN 2673 CALL GreenHyb_forFourier(this%Greens(iflavor1),Gomega=Gw(:,iflavor1)) 2674 ELSE 2675 CALL WARNALL("Ctqmc_getGreen : Gw is not valid ") 2676 CALL GreenHyb_forFourier(this%Greens(iflavor1),Wmax=this%Wmax) 2677 END IF 2678 ELSE 2679 CALL GreenHyb_forFourier(this%Greens(iflavor1),Wmax=this%Wmax) 2680 END IF 2681 END DO 2682 END SUBROUTINE Ctqmc_getGreen
ABINIT/m_Ctqmc/Ctqmc_getResult [ Functions ]
NAME
Ctqmc_getResult
FUNCTION
reduce everything to get the result of the simulation
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2105 SUBROUTINE Ctqmc_getResult(this) 2106 2107 2108 #ifdef HAVE_MPI1 2109 include 'mpif.h' 2110 #endif 2111 !Arguments ------------------------------------ 2112 TYPE(Ctqmc) , INTENT(INOUT) :: this 2113 !Local variables ------------------------------ 2114 INTEGER :: iflavor 2115 INTEGER :: flavors 2116 INTEGER :: itau 2117 INTEGER :: endDensity 2118 DOUBLE PRECISION :: inv_flavors 2119 DOUBLE PRECISION :: a 2120 DOUBLE PRECISION :: b 2121 DOUBLE PRECISION :: r 2122 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: alpha 2123 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: beta 2124 DOUBLE PRECISION, DIMENSION(1:2) :: TabX 2125 DOUBLE PRECISION, DIMENSION(1:2) :: TabY 2126 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: freqs 2127 INTEGER, ALLOCATABLE, DIMENSION(:) :: counts 2128 INTEGER, ALLOCATABLE, DIMENSION(:) :: displs 2129 INTEGER :: sp1 2130 INTEGER :: spAll 2131 INTEGER :: last 2132 INTEGER :: n1 2133 INTEGER :: n2 2134 INTEGER :: debut 2135 ! INTEGER :: fin 2136 #ifdef HAVE_MPI 2137 INTEGER :: ierr 2138 #endif 2139 DOUBLE PRECISION :: inv_size,sumh 2140 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: buffer 2141 TYPE(FFTHyb) :: FFTmrka 2142 2143 IF ( .NOT. this%done ) & 2144 CALL ERROR("Ctqmc_getResult : Simulation not run ") 2145 2146 flavors = this%flavors 2147 inv_flavors = 1.d0 / DBLE(flavors) 2148 2149 2150 inv_size = 1.d0 / DBLE(this%size) 2151 sp1 = 0 2152 spAll = 0 2153 2154 !#ifdef CTCtqmc_CHECK 2155 IF ( this%opt_check .GT. 0 ) THEN 2156 this%errorImpurity = ImpurityOperator_getError(this%Impurity) * inv_flavors 2157 this%errorBath = BathOperator_getError (this%Bath ) * inv_flavors 2158 END IF 2159 !#endif 2160 2161 MALLOC(alpha,(1,1)) 2162 MALLOC(beta,(1,1)) 2163 MALLOC(buffer,(1,1)) 2164 IF ( this%opt_noise .EQ. 1) THEN 2165 FREEIF(alpha) 2166 MALLOC(alpha,(1:this%samples+1,1:flavors)) 2167 FREEIF(beta) 2168 MALLOC(beta,(1:this%samples+1,1:flavors)) 2169 END IF 2170 2171 IF ( this%have_MPI .EQV. .TRUE.) THEN 2172 sp1 = this%samples+1 2173 spALL = sp1 + flavors + 6 2174 2175 !#ifdef CTCtqmc_ANALYSIS 2176 IF ( this%opt_analysis .EQ. 1 ) & 2177 spAll = spAll + 3*sp1 2178 IF ( this%opt_order .GT. 0 ) & 2179 spAll = spAll + this%opt_order 2180 IF ( this%opt_noise .EQ. 1 ) & 2181 spAll = spAll + 2*(this%samples + 1) 2182 !#endif 2183 2184 FREEIF(buffer) 2185 MALLOC(buffer,(1:spAll,1:MAX(2,flavors))) 2186 END IF 2187 2188 ! this%seg_added = this%seg_added * inv_flavors 2189 ! this%seg_removed = this%seg_removed * inv_flavors 2190 ! this%seg_sign = this%seg_sign * inv_flavors 2191 ! this%anti_added = this%anti_added * inv_flavors 2192 ! this%anti_removed = this%anti_removed * inv_flavors 2193 ! this%anti_sign = this%anti_sign * inv_flavors 2194 this%stats(:) = this%stats(:) * inv_flavors 2195 2196 DO iflavor = 1, flavors 2197 CALL GreenHyb_measHybrid(this%Greens(iflavor), this%Bath%M(iflavor), this%Impurity%Particles(iflavor), .TRUE.) 2198 CALL GreenHyb_getHybrid(this%Greens(iflavor)) 2199 ! Accumule les dernieres mesure de N 2200 this%measN(1,iflavor) = this%measN(1,iflavor) + this%measN(3,iflavor)*this%measN(4,iflavor) 2201 this%measN(2,iflavor) = this%measN(2,iflavor) + this%measN(4,iflavor) 2202 ! Reduction 2203 this%measN(1,iflavor) = this%measN(1,iflavor) / ( this%measN(2,iflavor) * this%beta ) 2204 ! Correction 2205 CALL GreenHyb_setN(this%Greens(iflavor), this%measN(1,iflavor)) 2206 !#ifdef CTCtqmc_ANALYSIS 2207 IF ( this%opt_order .GT. 0 ) & 2208 this%measPerturbation(: ,iflavor) = this%measPerturbation(:,iflavor) & 2209 / SUM(this%measPerturbation(:,iflavor)) 2210 IF ( this%opt_analysis .EQ. 1 ) THEN 2211 this%measCorrelation (:,1,iflavor) = this%measCorrelation (:,1,iflavor) & 2212 / SUM(this%measCorrelation (:,1,iflavor)) & 2213 * this%inv_dt 2214 this%measCorrelation (:,2,iflavor) = this%measCorrelation (:,2,iflavor) & 2215 / SUM(this%measCorrelation (:,2,iflavor)) & 2216 * this%inv_dt 2217 this%measCorrelation (:,3,iflavor) = this%measCorrelation (:,3,iflavor) & 2218 / SUM(this%measCorrelation (:,3,iflavor)) & 2219 * this%inv_dt 2220 END IF 2221 !#endif 2222 IF ( this%opt_noise .EQ. 1 ) THEN 2223 TabX(1) = DBLE(this%modNoise2) 2224 TabX(2) = DBLE(this%modNoise1) 2225 DO itau = 1, this%samples+1 2226 this%measNoiseG(itau,iflavor,2)%vec = -this%measNoiseG(itau,iflavor,2)%vec*this%inv_dt & 2227 /(this%beta*DBLE(this%modNoise2)) 2228 this%measNoiseG(itau,iflavor,1)%vec = -this%measNoiseG(itau,iflavor,1)%vec*this%inv_dt & 2229 /(this%beta*DBLE(this%modNoise1)) 2230 n2 = this%measNoiseG(itau,iflavor,2)%tail 2231 TabY(1) = Stat_deviation(this%measNoiseG(itau,iflavor,2)%vec(1:n2))!*SQRT(n2/(n2-1)) 2232 n1 = this%measNoiseG(itau,iflavor,1)%tail 2233 TabY(2) = Stat_deviation(this%measNoiseG(itau,iflavor,1)%vec(1:n1))!*SQRT(n1/(n1-1)) 2234 CALL Stat_powerReg(TabX,SQRT(2.d0*LOG(2.d0))*TabY,alpha(itau,iflavor),beta(itau,iflavor),r) 2235 ! ecart type -> 60% 2236 ! largeur a mi-hauteur d'une gaussienne -> sqrt(2*ln(2))*sigma 2237 END DO 2238 END IF 2239 2240 IF ( this%have_MPI .EQV. .TRUE. ) THEN 2241 buffer(1:sp1, iflavor) = this%Greens(iflavor)%oper(1:sp1) 2242 END IF 2243 END DO 2244 last = sp1 2245 2246 2247 this%measDE(:,:) = this%measDE(:,:) * DBLE(this%measurements) /(DBLE(this%sweeps)*this%beta) 2248 IF ( this%opt_histo .GT. 0 ) THEN 2249 this%occup_histo_time(:) = this%occup_histo_time(:) / INT(this%sweeps/this%measurements) 2250 END IF 2251 ! write(6,*) "=== Histogram of occupations for complete simulation ====",INT(this%sweeps/this%measurements) 2252 ! sumh=0 2253 ! do n1=1,this%flavors+1 2254 ! write(6,'(i4,f10.4)') n1-1, this%occup_histo_time(n1) 2255 ! sumh=sumh+this%occup_histo_time(n1) 2256 ! enddo 2257 ! write(6,*) "=================================",sumh 2258 2259 n1 = this%measNoise(1)%tail 2260 n2 = this%measNoise(2)%tail 2261 2262 ! On utilise freqs comme tableau de regroupement 2263 ! Gather de Noise1 2264 IF ( this%have_MPI .EQV. .TRUE. ) THEN 2265 MALLOC(counts,(1:this%size)) 2266 MALLOC(displs,(1:this%size)) 2267 FREEIF(freqs) 2268 MALLOC(freqs,(1:this%size*n1)) 2269 freqs = 0.d0 2270 freqs(n1*this%rank+1:n1*(this%rank+1)) = this%measNoise(1)%vec(1:n1) 2271 counts(:) = n1 2272 displs(:) = (/ ( iflavor*n1, iflavor=0, this%size-1 ) /) 2273 #ifdef HAVE_MPI 2274 CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_PRECISION, & 2275 freqs, counts, displs, & 2276 MPI_DOUBLE_PRECISION, this%MY_COMM, ierr) 2277 #endif 2278 n1 = this%size*n1 2279 CALL Vector_setSize(this%measNoise(1),n1) 2280 this%measNoise(1)%vec(1:n1) = freqs(:) 2281 ! Gather de Noise2 2282 FREE(freqs) 2283 MALLOC(freqs,(1:this%size*n2)) 2284 freqs = 0.d0 2285 freqs(n2*this%rank+1:n2*(this%rank+1)) = this%measNoise(2)%vec(1:n2) 2286 counts(:) = n2 2287 displs(:) = (/ ( iflavor*n2, iflavor=0, this%size-1 ) /) 2288 #ifdef HAVE_MPI 2289 CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DATATYPE_NULL, & 2290 freqs, counts, displs, & 2291 MPI_DOUBLE_PRECISION, this%MY_COMM, ierr) 2292 #endif 2293 n2 = this%size*n2 2294 CALL Vector_setSize(this%measNoise(2),n2) 2295 this%measNoise(2)%vec(1:n2) = freqs(:) 2296 FREE(counts) 2297 FREE(displs) 2298 FREE(freqs) 2299 END IF 2300 !n1 = this%measNoise(1)%tail 2301 !n2 = this%measNoise(2)%tail 2302 2303 ! Transformation des paquets pour que ca fit a CTQMC_SLICE(1|2) 2304 IF ( n1 .GT. CTQMC_SLICE1 ) THEN 2305 itau = n1/CTQMC_SLICE1 2306 MALLOC(freqs,(1:n1/itau)) 2307 DO debut=1, n1/itau 2308 freqs(debut)=SUM(this%measNoise(1)%vec((debut-1)*itau+1:itau*debut)) 2309 END DO 2310 freqs(:) = freqs(:)/DBLE(itau) 2311 this%modNoise1 = this%modNoise1*itau 2312 n1 = n1/itau 2313 CALL Vector_setSize(this%measNoise(1),n1) 2314 this%measNoise(1)%vec(1:n1) = freqs(:) 2315 FREE(freqs) 2316 END IF 2317 IF ( n2 .GT. CTQMC_SLICE1*CTQMC_SLICE2 ) THEN 2318 itau = n2/(CTQMC_SLICE1*CTQMC_SLICE2) 2319 MALLOC(freqs,(1:n2/itau)) 2320 DO debut=1, n2/itau 2321 freqs(debut)=SUM(this%measNoise(2)%vec((debut-1)*itau+1:itau*debut)) 2322 END DO 2323 freqs(:) = freqs(:)/DBLE(itau) 2324 this%modNoise2 = this%modNoise2*itau 2325 n2 = n2/itau 2326 CALL Vector_setSize(this%measNoise(2),n2) 2327 this%measNoise(2)%vec(1:n2) = freqs(:) 2328 FREE(freqs) 2329 END IF 2330 ! On peut s'amuser avec nos valeur d'energies 2331 !MALLOC(TabX,(1:20)) 2332 !MALLOC(TabY,(1:20)) 2333 2334 TabX(1) = DBLE(this%modNoise2) 2335 TabX(2) = DBLE(this%modNoise1) 2336 2337 ! Il faut calculer pour chaque modulo 10 ecarts type sur les donnes acquises 2338 this%measNoise(1)%vec(1:n1) = this%measNoise(1)%vec(1:n1)/(this%beta*DBLE(this%modNoise1))*DBLE(this%measurements) 2339 this%measNoise(2)%vec(1:n2) = this%measNoise(2)%vec(1:n2)/(this%beta*DBLE(this%modNoise2))*DBLE(this%measurements) 2340 ! CALL Vector_print(this%measNoise(1),this%rank+70) 2341 ! CALL Vector_print(this%measNoise(2),this%rank+50) 2342 ! DO iflavor=1,10 2343 ! debut = (iflavor-1)*n2/10+1 2344 ! fin = iflavor*n2/10 2345 ! TabY(iflavor) = Stat_deviation(this%measNoise(2)%vec(debut:fin)) 2346 ! debut = (iflavor-1)*n1/10+1 2347 ! fin = iflavor*n1/10 2348 ! TabY(10+iflavor) = Stat_deviation(this%measNoise(1)%vec(debut:fin)) 2349 ! END DO 2350 !! TabY(1:n) = (this%measNoise(2)%vec(1:n) & 2351 !! ) 2352 !! !/(this%beta*DBLE(this%modNoise2))*DBLE(this%measurements) & 2353 !! !- this%measDE(1,1)) 2354 !! TabY(this%measNoise(2)%tail+1:n+this%measNoise(2)%tail) = (this%measNoise(1)%vec(1:n) & 2355 !! ) 2356 !! ! /(this%beta*DBLE(this%modNoise1))*DBLE(this%measurements) & 2357 !! ! - this%measDE(1,1)) 2358 ! IF ( this%rank .EQ. 0 ) THEN 2359 ! DO iflavor=1,20 2360 ! write(45,*) TabX(iflavor), TabY(iflavor) 2361 ! END DO 2362 ! END IF 2363 ! 2364 2365 2366 TabY(1) = Stat_deviation(this%measNoise(2)%vec(1:n2))!*SQRT(n2/(n2-1)) 2367 !! write(this%rank+10,*) TabX(2) 2368 !! write(this%rank+40,*) TabX(1) 2369 !! CALL Vector_print(this%measNoise(1),this%rank+10) 2370 !! CALL Vector_print(this%measNoise(2),this%rank+40) 2371 !! CLOSE(this%rank+10) 2372 !! CLOSE(this%rank+40) 2373 TabY(2) = Stat_deviation(this%measNoise(1)%vec(1:n1))!*SQRT(n1/(n1-1)) 2374 !! ! Ecart carre moyen ~ ecart type mais non biaise. Serait moins precis. Aucun 2375 ! impact sur la pente, juste sur l'ordonnee a l'origine. 2376 2377 CALL Stat_powerReg(TabX,SQRT(2.d0*LOG(2.d0))*TabY,a,b,r) 2378 ! FREE(TabX) 2379 ! FREE(TabY) 2380 ! ecart type -> 60% 2381 ! largeur a mi-hauteur d'une gaussienne -> sqrt(2*ln(2))*sigma 2382 2383 !this%measDE(1,1) = SUM(this%measNoise(1)%vec(1:this%measNoise(1)%tail))/(DBLE(this%measNoise(1)%tail*this%modNoise1)*this%beta) 2384 !this%measDE(2:flavors,1:flavors) = this%measDE(2:flavors,1:flavors) /(DBLE(this%sweeps)*this%beta) 2385 CALL ImpurityOperator_getErrorOverlap(this%Impurity,this%measDE) 2386 ! Add the difference between true calculation and quick calculation of the 2387 ! last sweep overlap to measDE(2,2) 2388 !this%measDE = this%measDE * DBLE(this%measurements) 2389 IF ( this%have_MPI .EQV. .TRUE. ) THEN 2390 IF ( this%opt_analysis .EQ. 1 ) THEN 2391 buffer(last+1:last+sp1,:) = this%measCorrelation(:,1,:) 2392 last = last + sp1 2393 buffer(last+1:last+sp1,:) = this%measCorrelation(:,2,:) 2394 last = last + sp1 2395 buffer(last+1:last+sp1,:) = this%measCorrelation(:,3,:) 2396 last = last + sp1 2397 END IF 2398 IF ( this%opt_order .GT. 0 ) THEN 2399 buffer(last+1:last+this%opt_order, :) = this%measPerturbation(:,:) 2400 last = last + this%opt_order 2401 END IF 2402 IF ( this%opt_noise .EQ. 1 ) THEN 2403 buffer(last+1:last+this%samples+1,:) = alpha(:,:) 2404 last = last + this%samples + 1 2405 buffer(last+1:last+this%samples+1,:) = beta(:,:) 2406 last = last + this%samples + 1 2407 END IF 2408 ! this%measDE(2,2) = a*EXP(b*LOG(DBLE(this%sweeps*this%size))) 2409 buffer(spall-(flavors+5):spAll-6,:) = this%measDE(:,:) 2410 ! buffer(spAll ,1) = this%seg_added 2411 ! buffer(spAll-1,1) = this%seg_removed 2412 ! buffer(spAll-2,1) = this%seg_sign 2413 ! buffer(spAll ,2) = this%anti_added 2414 ! buffer(spAll-1,2) = this%anti_removed 2415 ! buffer(spAll-2,2) = this%anti_sign 2416 buffer(spAll ,1) = this%stats(1) 2417 buffer(spAll-1,1) = this%stats(2) 2418 buffer(spAll-2,1) = this%stats(3) 2419 buffer(spAll ,2) = this%stats(4) 2420 buffer(spAll-1,2) = this%stats(5) 2421 buffer(spAll-2,2) = this%stats(6) 2422 buffer(spAll-3,1) = this%swap 2423 buffer(spAll-3,2) = DBLE(this%modGlobalMove(2)) 2424 buffer(spAll-4,1) = a 2425 buffer(spAll-4,2) = b 2426 !#ifdef CTCtqmc_CHECK 2427 buffer(spAll-5,1) = this%errorImpurity 2428 buffer(spAll-5,2) = this%errorBath 2429 !#endif 2430 2431 #ifdef HAVE_MPI 2432 CALL MPI_ALLREDUCE(MPI_IN_PLACE, buffer, spAll*flavors, & 2433 MPI_DOUBLE_PRECISION, MPI_SUM, this%MY_COMM, ierr) 2434 CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%runTime, 1, MPI_DOUBLE_PRECISION, MPI_MAX, & 2435 this%MY_COMM, ierr) 2436 IF ( this%opt_histo .GT. 0 ) THEN 2437 CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%occup_histo_time, flavors+1, MPI_DOUBLE_PRECISION, MPI_SUM, & 2438 this%MY_COMM, ierr) 2439 END IF 2440 #endif 2441 2442 2443 buffer = buffer * inv_size 2444 this%measDE(:,:) = buffer(spall-(flavors+5):spAll-6,:) 2445 ! this%seg_added = buffer(spAll ,1) 2446 ! this%seg_removed = buffer(spAll-1,1) 2447 ! this%seg_sign = buffer(spAll-2,1) 2448 ! this%anti_added = buffer(spAll ,2) 2449 ! this%anti_removed = buffer(spAll-1,2) 2450 ! this%anti_sign = buffer(spAll-2,2) 2451 this%stats(1) = buffer(spAll ,1) 2452 this%stats(2) = buffer(spAll-1,1) 2453 this%stats(3) = buffer(spAll-2,1) 2454 this%stats(4) = buffer(spAll ,2) 2455 this%stats(5) = buffer(spAll-1,2) 2456 this%stats(6) = buffer(spAll-2,2) 2457 this%swap = buffer(spAll-3,1) 2458 this%modGlobalMove(2) = NINT(buffer(spAll-3,2)) 2459 a = buffer(spAll-4,1) 2460 b = buffer(spAll-4,2) 2461 !!#ifdef CTCtqmc_CHECK 2462 this%errorImpurity= buffer(spAll-5,1) 2463 this%errorBath = buffer(spAll-5,2) 2464 !#endif 2465 2466 DO iflavor = 1, flavors 2467 this%Greens(iflavor)%oper = buffer(1:sp1 , iflavor) 2468 END DO 2469 last = sp1 2470 IF ( this%opt_analysis .EQ. 1 ) THEN 2471 this%measCorrelation(:,1,:) = buffer(last+1:last+sp1,:) 2472 last = last + sp1 2473 this%measCorrelation(:,2,:) = buffer(last+1:last+sp1,:) 2474 last = last + sp1 2475 this%measCorrelation(:,3,:) = buffer(last+1:last+sp1,:) 2476 last = last + sp1 2477 END IF 2478 IF ( this%opt_order .GT. 0 ) THEN 2479 this%measPerturbation(:,:) = buffer(last+1:last+this%opt_order, :) 2480 last = last + this%opt_order 2481 END IF 2482 IF ( this%opt_noise .EQ. 1 ) THEN 2483 alpha(:,:) = buffer(last+1:last+this%samples+1,:) 2484 last = last + this%samples + 1 2485 beta(:,:) = buffer(last+1:last+this%samples+1,:) 2486 last = last + this%samples + 1 2487 END IF 2488 END IF 2489 DO iflavor = 1, flavors 2490 ! complete DE this 2491 this%measDE(iflavor, iflavor+1:flavors) = this%measDE(iflavor+1:flavors,iflavor) 2492 END DO 2493 FREE(buffer) 2494 2495 IF ( this%opt_spectra .GE. 1 ) THEN 2496 endDensity = SIZE(this%density,2) 2497 IF ( this%density(1,endDensity) .EQ. -1.d0 ) & 2498 endDensity = endDensity - 1 2499 CALL FFTHyb_init(FFTmrka,endDensity,DBLE(this%thermalization)/DBLE(this%measurements*this%opt_spectra)) 2500 ! Not very Beauty 2501 MALLOC(freqs,(1:FFTmrka%size/2)) 2502 DO iflavor = 1, flavors 2503 ! mean value is removed to supress the continue composent 2504 CALL FFTHyb_setData(FFTmrka,this%density(iflavor,1:endDensity)/this%beta+this%Greens(iflavor)%oper(this%samples+1)) 2505 CALL FFTHyb_run(FFTmrka,1) 2506 CALL FFTHyb_getData(FFTmrka,endDensity,this%density(iflavor,:),freqs) 2507 END DO 2508 this%density(flavors+1,:) = -1.d0 2509 this%density(flavors+1,1:FFTmrka%size/2) = freqs 2510 CALL FFTHyb_destroy(FFTmrka) 2511 FREE(freqs) 2512 END IF 2513 2514 this%a_Noise = a 2515 this%b_Noise = b 2516 IF ( this%opt_noise .EQ. 1 ) THEN 2517 this%abNoiseG(1,:,:) = alpha 2518 this%abNoiseG(2,:,:) = beta 2519 END IF 2520 FREE(alpha) 2521 FREE(beta) 2522 2523 IF ( this%opt_histo .GT. 0 ) THEN 2524 write(this%ostream,*) "=== Histogram of occupations for complete simulation ====" 2525 ! write(6,*) "sumh over procs", sumh 2526 sumh=0 2527 do n1=1,this%flavors+1 2528 write(this%ostream,'(i4,f10.4)') n1-1, this%occup_histo_time(n1)/float(this%size) 2529 sumh=sumh+this%occup_histo_time(n1)/float(this%size) 2530 enddo 2531 write(this%ostream,'(a,f10.4)') " all" , sumh 2532 write(this%ostream,*) "=================================" 2533 ENDIF 2534 2535 END SUBROUTINE Ctqmc_getResult
ABINIT/m_Ctqmc/Ctqmc_init [ Functions ]
NAME
Ctqmc_init
FUNCTION
Initialize the type Ctqmc Allocate all the non optional variables
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc ostream=where to write istream=where to read the input parameters if so bFile=do we read in istream ? MY_COMM=mpi communicator for the CTQMC iBuffer=input parameters if bFile is false
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
306 SUBROUTINE Ctqmc_init(this, ostream, istream, bFile, MY_COMM, iBuffer) 307 308 309 #ifdef HAVE_MPI1 310 include 'mpif.h' 311 #endif 312 !Arguments ------------------------------------ 313 TYPE(Ctqmc), INTENT(INOUT) :: this 314 INTEGER , INTENT(IN ) :: ostream 315 INTEGER , INTENT(IN ) :: istream 316 LOGICAL , INTENT(IN ) :: bFile 317 DOUBLE PRECISION, DIMENSION(1:9), OPTIONAL, INTENT(IN) :: iBuffer 318 INTEGER , OPTIONAL, INTENT(IN ) :: MY_COMM 319 !Local variables ------------------------------ 320 #ifdef HAVE_MPI 321 INTEGER :: ierr 322 #endif 323 INTEGER :: iflavor 324 #ifdef __GFORTRAN__ 325 ! INTEGER :: pid 326 ! CHARACTER(LEN=5) :: Cpid 327 ! 328 #endif 329 DOUBLE PRECISION, DIMENSION(1:9) :: buffer 330 331 this%ostream = ostream 332 this%istream = istream 333 334 ! --- RENICE --- 335 !#ifdef __GFORTRAN__ 336 ! pid = GetPid() 337 ! WRITE(Cpid,'(I5)') pid 338 ! CALL SYSTEM('renice +19 '//TRIM(ADJUSTL(Cpid))//' > /dev/null') 339 !#endif 340 !! --- RENICE --- 341 342 IF ( PRESENT(MY_COMM)) THEN 343 #ifdef HAVE_MPI 344 this%have_MPI = .TRUE. 345 this%MY_COMM = MY_COMM 346 CALL MPI_Comm_rank(this%MY_COMM, this%rank, ierr) 347 CALL MPI_Comm_size(this%MY_COMM, this%size, ierr) 348 #else 349 CALL WARN("Ctqmc_init : MPI is not used ") 350 this%have_MPI = .FALSE. 351 this%MY_COMM = -1 352 this%rank = 0 353 this%size = 1 354 #endif 355 ELSE 356 this%have_MPI = .FALSE. 357 this%MY_COMM = -1 358 this%rank = 0 359 this%size = 1 360 END IF 361 362 !IF ( this%rank .EQ. 0 ) THEN 363 ! WRITE(ostream,'(A20)') 'Job reniced with +19' 364 !CALL FLUSH(ostream) 365 !END IF 366 367 IF ( bFile .EQV. .TRUE. ) THEN 368 IF ( this%rank .EQ. 0 ) THEN 369 370 READ(istream,*) buffer(1) !iseed 371 READ(istream,*) buffer(2) !this%sweeps 372 READ(istream,*) buffer(3) !this%thermalization 373 READ(istream,*) buffer(4) !this%measurements 374 READ(istream,*) buffer(5) !this%flavors 375 READ(istream,*) buffer(6) !this%samples 376 READ(istream,*) buffer(7) !this%beta 377 READ(istream,*) buffer(8) !U 378 READ(istream,*) buffer(9) !iTech 379 !READ(istream,*) buffer(9) !Wmax 380 !#ifdef CTCtqmc_ANALYSIS 381 !READ(istream,*) buffer(10) !order 382 !#endif 383 END IF 384 385 #ifdef HAVE_MPI 386 IF ( this%have_MPI .EQV. .TRUE. ) & 387 CALL MPI_Bcast(buffer, 9, MPI_DOUBLE_PRECISION, 0, & 388 this%MY_COMM, ierr) 389 #endif 390 ELSE IF ( PRESENT(iBuffer) ) THEN 391 buffer(1:9) = iBuffer(1:9) 392 ELSE 393 CALL ERROR("Ctqmc_init : No input parameters ") 394 END IF 395 396 CALL Ctqmc_setParameters(this, buffer) 397 398 CALL Ctqmc_allocateAll(this) 399 400 DO iflavor = 1, this%flavors 401 CALL GreenHyb_init(this%Greens(iflavor),this%samples, this%beta, iTech=INT(buffer(9)),MY_COMM=this%MY_COMM) 402 END DO 403 404 405 ! this%seg_added = 0.d0 406 ! this%anti_added = 0.d0 407 ! this%seg_removed = 0.d0 408 ! this%anti_removed = 0.d0 409 ! this%seg_sign = 0.d0 410 ! this%anti_sign = 0.d0 411 this%stats(:) = 0.d0 412 this%swap = 0.d0 413 this%runTime = 0.d0 414 415 CALL Vector_init(this%measNoise(1),this%sweeps/this%modNoise1) 416 CALL Vector_init(this%measNoise(2),(this%sweeps/this%modNoise1+1)*CTQMC_SLICE2) 417 !CALL Vector_init(this%measNoise(3),101) 418 !CALL Vector_init(this%measNoise(4),101) 419 420 this%set = this%para .AND. this%inF 421 this%done = .FALSE. 422 this%init = .TRUE. 423 424 !#ifdef CTCtqmc_CHECK 425 this%errorImpurity = 0.d0 426 this%errorBath = 0.d0 427 !#endif 428 END SUBROUTINE Ctqmc_init
ABINIT/m_Ctqmc/Ctqmc_loop [ Functions ]
NAME
Ctqmc_loop
FUNCTION
Definition the main loop of the CT-QMC
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc itotal=number of sweeps to perform : thermalization or sweeps ilatex=unit of file to write movie if so
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1464 SUBROUTINE Ctqmc_loop(this,itotal,ilatex) 1465 1466 !Arguments ------------------------------------ 1467 TYPE(Ctqmc), INTENT(INOUT) :: this 1468 INTEGER , INTENT(IN ) :: itotal 1469 INTEGER , INTENT(IN ) :: ilatex 1470 !Local variables ------------------------------ 1471 LOGICAL :: updated 1472 LOGICAL :: updated_seg 1473 LOGICAL, DIMENSION(:), ALLOCATABLE :: updated_swap 1474 1475 INTEGER :: flavors 1476 INTEGER :: measurements 1477 INTEGER :: modNoise1 1478 INTEGER :: modNoise2 1479 INTEGER :: modGlobalMove 1480 INTEGER :: sp1 1481 INTEGER :: itau 1482 INTEGER :: ind 1483 INTEGER :: endDensity 1484 INTEGER :: indDensity 1485 INTEGER :: swapUpdate1 1486 INTEGER :: swapUpdate2 1487 INTEGER :: old_percent 1488 INTEGER :: new_percent 1489 INTEGER :: ipercent 1490 INTEGER :: iflavor 1491 INTEGER :: isweep 1492 1493 DOUBLE PRECISION :: cpu_time1 1494 DOUBLE PRECISION :: cpu_time2 1495 DOUBLE PRECISION :: NRJ_old1 1496 DOUBLE PRECISION :: NRJ_old2 1497 DOUBLE PRECISION :: NRJ_new 1498 DOUBLE PRECISION :: total 1499 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_old1 1500 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_old2 1501 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_new 1502 1503 CALL CPU_TIME(cpu_time1) 1504 1505 flavors = this%flavors 1506 measurements = this%measurements 1507 modNoise1 = this%modNoise1 1508 modNoise2 = this%modNoise2 1509 modGlobalMove = this%modGlobalMove(1) 1510 sp1 = this%samples+1 1511 IF ( this%opt_histo .GT. 0 ) THEN 1512 this%occup_histo_time= 0.d0 1513 END IF 1514 old_percent = 0 1515 1516 MALLOC(updated_swap,(1:flavors)) 1517 updated_swap(:) = .FALSE. 1518 1519 NRJ_old1 = 0.d0 1520 NRJ_old2 = 0.d0 1521 NRJ_new = 0.d0 1522 1523 MALLOC(gtmp_new,(1,1)) 1524 gtmp_new = 0.d0 1525 MALLOC(gtmp_old1,(1,1)) 1526 gtmp_old1 = 0.d0 1527 MALLOC(gtmp_old2,(1,1)) 1528 gtmp_old2 = 0.d0 1529 1530 endDensity = SIZE(this%density,2) 1531 1532 IF ( this%opt_noise .GT. 0 ) THEN 1533 FREEIF(gtmp_new) 1534 MALLOC(gtmp_new,(1:sp1,1:flavors)) 1535 FREEIF(gtmp_old1) 1536 MALLOC(gtmp_old1,(1:sp1,1:flavors)) 1537 FREEIF(gtmp_old2) 1538 MALLOC(gtmp_old2,(1:sp1,1:flavors)) 1539 END IF 1540 1541 IF ( this%rank .EQ. 0 ) THEN 1542 WRITE(this%ostream, '(1x,103A)') & 1543 "|----------------------------------------------------------------------------------------------------|" 1544 WRITE(this%ostream,'(1x,A)', ADVANCE="NO") "|" 1545 END IF 1546 1547 total = DBLE(itotal) 1548 1549 indDensity = 1 1550 DO isweep = 1, itotal 1551 DO iflavor = 1, flavors 1552 ImpurityOperator_QuickActivation(this%Impurity,iflavor) 1553 BathOperator_QuickActivation(this%Bath,iflavor) 1554 1555 CALL Ctqmc_tryAddRemove(this,updated_seg) 1556 updated = updated_seg .OR. updated_swap(iflavor) 1557 updated_swap(iflavor) = .FALSE. 1558 1559 CALL GreenHyb_measHybrid(this%Greens(iflavor), this%Bath%M(iflavor), this%Impurity%Particles(iflavor), updated) 1560 CALL Ctqmc_measN (this, iflavor, updated) 1561 IF ( this%opt_analysis .EQ. 1 ) & 1562 CALL Ctqmc_measCorrelation (this, iflavor) 1563 IF ( this%opt_order .GT. 0 ) & 1564 CALL Ctqmc_measPerturbation(this, iflavor) 1565 END DO 1566 1567 IF ( MOD(isweep,modGlobalMove) .EQ. 0 ) THEN 1568 CALL Ctqmc_trySwap(this,swapUpdate1, swapUpdate2) 1569 IF ( swapUpdate1 .NE. 0 .AND. swapUpdate2 .NE. 0 ) THEN 1570 updated_swap(swapUpdate1) = .TRUE. 1571 updated_swap(swapUpdate2) = .TRUE. 1572 END IF 1573 END IF 1574 1575 IF ( MOD(isweep,measurements) .EQ. 0 ) THEN 1576 CALL ImpurityOperator_measDE(this%Impurity,this%measDE) 1577 IF ( this%opt_spectra .GE. 1) THEN 1578 IF ( measurements*this%opt_spectra /= 0) THEN 1579 IF ( MOD(isweep,measurements*this%opt_spectra) .EQ. 0 ) THEN 1580 this%density(1:flavors,indDensity) = this%measN(3,1:flavors) 1581 indDensity = indDensity+1 1582 END IF 1583 END IF 1584 END IF 1585 END IF 1586 1587 IF ( MOD(isweep,measurements) .EQ. 0 ) THEN 1588 IF ( this%opt_histo .GT. 0 ) THEN 1589 CALL ImpurityOperator_occup_histo_time(this%Impurity,this%occup_histo_time) 1590 ENDIF 1591 ENDIF 1592 1593 IF ( MOD(isweep, modNoise1) .EQ. 0 ) THEN 1594 !modNext = isweep + modNoise2 1595 NRJ_new = (SUM(this%measDE(:,:))-this%measDE(1,1))*0.5d0 ! double occupation, avoid stat with 0 for U=J=0 1596 CALL Vector_pushBack(this%measNoise(1),NRJ_new - NRJ_old1) 1597 NRJ_old1 = NRJ_new 1598 1599 !! Try to limit accumulation error 1600 CALL ImpurityOperator_cleanOverlaps(this%Impurity) 1601 1602 IF ( this%opt_noise .EQ. 1 ) THEN 1603 DO iflavor = 1, flavors 1604 DO ind = 1, this%Greens(iflavor)%this%tail 1605 itau = this%Greens(iflavor)%this%listINT(ind) 1606 gtmp_new(itau,iflavor) = this%Greens(iflavor)%oper(itau) & 1607 +this%Greens(iflavor)%this%listDBLE(ind)*DBLE(this%Greens(iflavor)%factor) 1608 END DO 1609 DO itau = 1, sp1 1610 CALL Vector_pushBack(this%measNoiseG(itau,iflavor,1), gtmp_new(itau,iflavor) - gtmp_old1(itau,iflavor)) 1611 gtmp_old1(itau,iflavor) = gtmp_new(itau,iflavor) 1612 END DO 1613 END DO 1614 END IF 1615 END IF 1616 1617 IF ( MOD(isweep,modNoise2) .EQ. 0 ) THEN 1618 NRJ_new = (SUM(this%measDE(:,:))-this%measDE(1,1))*0.5d0 ! double occupation, avoid stat with 0 for U=J=0 1619 CALL Vector_pushBack(this%measNoise(2),NRJ_new - NRJ_old2) 1620 NRJ_old2 = NRJ_new 1621 IF ( this%opt_noise .EQ. 1 ) THEN 1622 DO iflavor = 1, flavors 1623 DO ind = 1, this%Greens(iflavor)%this%tail 1624 itau = this%Greens(iflavor)%this%listINT(ind) 1625 gtmp_new(itau,iflavor) = this%Greens(iflavor)%oper(itau) & 1626 +this%Greens(iflavor)%this%listDBLE(ind)*this%Greens(iflavor)%factor 1627 END DO 1628 DO itau = 1, sp1 1629 CALL Vector_pushBack(this%measNoiseG(itau,iflavor,2), gtmp_new(itau,iflavor) - gtmp_old2(itau,iflavor)) 1630 gtmp_old2(itau,iflavor) = gtmp_new(itau,iflavor) 1631 END DO 1632 END DO 1633 END IF 1634 1635 IF ( this%rank .EQ. 0 ) THEN 1636 new_percent = CEILING(DBLE(isweep)*100.d0/DBLE(itotal)) 1637 DO ipercent = old_percent+1, new_percent 1638 WRITE(this%ostream,'(A)',ADVANCE="NO") "-" 1639 END DO 1640 old_percent = new_percent 1641 END IF 1642 END IF 1643 1644 IF ( this%opt_movie .EQ. 1 ) THEN 1645 WRITE(ilatex,'(A11,I9)') "%iteration ", isweep 1646 CALL ImpurityOperator_printLatex(this%Impurity,ilatex,isweep) 1647 END IF 1648 1649 END DO 1650 1651 IF ( this%rank .EQ. 0 ) THEN 1652 DO ipercent = old_percent+1, 100 1653 WRITE(this%ostream,'(A)',ADVANCE="NO") "-" 1654 END DO 1655 WRITE(this%ostream,'(A)') "|" 1656 END IF 1657 1658 FREE(gtmp_new) 1659 FREE(gtmp_old1) 1660 FREE(gtmp_old2) 1661 FREE(updated_swap) 1662 1663 IF ( this%opt_spectra .GE. 1 .AND. itotal .EQ. this%sweeps ) THEN 1664 IF ( endDensity .NE. indDensity-1 ) THEN 1665 this%density(:,endDensity) = -1.d0 1666 END IF 1667 END IF 1668 1669 CALL CPU_TIME(cpu_time2) 1670 1671 this%runTime = (cpu_time2 - cpu_time1)*1.05d0 ! facteur arbitraire de correction 1672 END SUBROUTINE Ctqmc_loop
ABINIT/m_Ctqmc/Ctqmc_measCorrelation [ Functions ]
NAME
Ctqmc_measCorrelation
FUNCTION
measure all correlations in times for a flavor
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc iflavor=the flavor to measure
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1977 SUBROUTINE Ctqmc_measCorrelation(this, iflavor) 1978 1979 !Arguments ------------------------------------ 1980 TYPE(Ctqmc) , INTENT(INOUT) :: this 1981 !TYPE(ImpurityOperator), INTENT(IN ) :: impurity 1982 INTEGER , INTENT(IN ) :: iflavor 1983 !Local variables ------------------------------ 1984 INTEGER :: iCdag 1985 INTEGER :: iCdagBeta 1986 INTEGER :: iC 1987 INTEGER :: index 1988 INTEGER :: size 1989 DOUBLE PRECISION :: tC 1990 DOUBLE PRECISION :: tCdag 1991 !DOUBLE PRECISION :: time 1992 DOUBLE PRECISION :: inv_dt 1993 DOUBLE PRECISION :: beta 1994 1995 IF ( .NOT. this%set ) & 1996 CALL ERROR("Ctqmc_measCorrelation : QMC not set ") 1997 1998 size = this%impurity%particles(this%impurity%activeFlavor)%tail 1999 beta = this%beta 2000 2001 IF ( size .EQ. 0 ) RETURN 2002 2003 inv_dt = this%inv_dt 2004 2005 DO iCdag = 1, size ! first segments 2006 tCdag = this%impurity%particles(this%impurity%activeFlavor)%list(iCdag,Cdag_) 2007 tC = this%impurity%particles(this%impurity%activeFlavor)%list(iCdag,C_ ) 2008 index = INT( ( (tC - tCdag) * inv_dt ) + .5d0 ) + 1 2009 this%measCorrelation(index,1,iflavor) = this%measCorrelation(index,1,iflavor) + 1.d0 2010 MODCYCLE(iCdag+1,size,iCdagBeta) 2011 index = INT( ( ( & 2012 this%impurity%particles(this%impurity%activeFlavor)%list(iCdagBeta,Cdag_) - tC & 2013 + AINT(DBLE(iCdag)/DBLE(size))*beta & 2014 ) * inv_dt ) + .5d0 ) + 1 2015 IF ( index .LT. 1 .OR. index .GT. this%samples+1 ) THEN 2016 CALL WARN("Ctqmc_measCorrelation : bad index line 1095 ") 2017 ELSE 2018 this%measCorrelation(index,2,iflavor) = this%measCorrelation(index,2,iflavor) + 1.d0 2019 END IF 2020 ! DO iC = 1, size 2021 ! tC = impurity%particles(impurity%activeFlavor)%list(C_,iC) 2022 ! time = tC - tCdag 2023 ! IF ( time .LT. 0.d0 ) time = time + beta 2024 ! index = INT( ( time * inv_dt ) + .5d0 ) + 1 2025 ! this%measCorrelation(index,3,iflavor) = this%measCorrelation(index,3,iflavor) + 1.d0 2026 ! END DO 2027 DO iC = 1, size! this%Greens(iflavor)%index_old%tail 2028 this%measCorrelation(this%Greens(iflavor)%this%listINT(iC+(iCdag-1)*size),3,iflavor) = & 2029 this%measCorrelation(this%Greens(iflavor)%this%listINT(iC+(iCdag-1)*size),3,iflavor) + 1.d0 2030 END DO 2031 END DO 2032 2033 END SUBROUTINE Ctqmc_measCorrelation
ABINIT/m_Ctqmc/Ctqmc_measN [ Functions ]
NAME
Ctqmc_measN
FUNCTION
measure the number of electron
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc iflavor=which flavor to measure updated=something has changed since last time
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1928 SUBROUTINE Ctqmc_measN(this, iflavor, updated) 1929 1930 !Arguments ------------------------------------ 1931 TYPE(Ctqmc) , INTENT(INOUT) :: this 1932 !TYPE(ImpurityOperator), INTENT(IN ) :: impurity 1933 INTEGER , INTENT(IN ) :: iflavor 1934 LOGICAL , INTENT(IN ) :: updated 1935 1936 ! IF ( .NOT. this%set ) & 1937 ! CALL ERROR("Ctqmc_measN : QMC not set ") 1938 1939 1940 IF ( updated .EQV. .TRUE. ) THEN 1941 this%measN(1,iflavor) = this%measN(1,iflavor) + this%measN(3,iflavor)*this%measN(4,iflavor) 1942 this%measN(2,iflavor) = this%measN(2,iflavor) + this%measN(4,iflavor) 1943 this%measN(3,iflavor) = ImpurityOperator_measN(this%impurity) 1944 this%measN(4,iflavor) = 1.d0 1945 ELSE 1946 this%measN(4,iflavor) = this%measN(4,iflavor) + 1.d0 1947 END IF 1948 END SUBROUTINE Ctqmc_measN
ABINIT/m_Ctqmc/Ctqmc_measPerturbation [ Functions ]
NAME
Ctqmc_measPerturbation
FUNCTION
measure perturbation order
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc iflavor=the flavor to measure
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2061 SUBROUTINE Ctqmc_measPerturbation(this, iflavor) 2062 2063 !Arguments ------------------------------------ 2064 TYPE(Ctqmc) , INTENT(INOUT) :: this 2065 !TYPE(ImpurityOperator), INTENT(IN ) :: impurity 2066 INTEGER , INTENT(IN ) :: iflavor 2067 !Local variables ------------------------------ 2068 INTEGER :: index 2069 2070 IF ( .NOT. this%set ) & 2071 CALL ERROR("Ctqmc_measiPerturbation : QMC not set ") 2072 2073 index = this%impurity%particles(this%impurity%activeFlavor)%tail + 1 2074 IF ( index .LE. this%opt_order ) & 2075 this%measPerturbation(index,iflavor) = this%measPerturbation(index,iflavor) + 1.d0 2076 2077 END SUBROUTINE Ctqmc_measPerturbation
ABINIT/m_Ctqmc/Ctqmc_printAll [ Functions ]
NAME
Ctqmc_printAll
FUNCTION
print different functions computed during the simulation
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2799 SUBROUTINE Ctqmc_printAll(this) 2800 2801 !Arguments ------------------------------------ 2802 TYPE(Ctqmc), INTENT(INOUT) :: this 2803 2804 IF ( .NOT. this%done ) & 2805 CALL WARNALL("Ctqmc_printAll : Simulation not run ") 2806 2807 CALL Ctqmc_printQMC(this) 2808 2809 CALL Ctqmc_printGreen(this) 2810 2811 CALL Ctqmc_printD(this) 2812 2813 ! CALL Ctqmc_printE(this) 2814 2815 !#ifdef CTCtqmc_ANALYSIS 2816 CALL Ctqmc_printPerturbation(this) 2817 2818 CALL Ctqmc_printCorrelation(this) 2819 !#endif 2820 2821 CALL Ctqmc_printSpectra(this) 2822 2823 END SUBROUTINE Ctqmc_printAll
ABINIT/m_Ctqmc/Ctqmc_printCorrelation [ Functions ]
NAME
Ctqmc_printCorrelation
FUNCTION
print correlation fonctions
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3207 SUBROUTINE Ctqmc_printCorrelation(this, oFileIn) 3208 3209 !Arguments ------------------------------------ 3210 TYPE(Ctqmc) , INTENT(IN) :: this 3211 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3212 !Local variables ------------------------------ 3213 INTEGER :: oFile 3214 INTEGER :: itime 3215 INTEGER :: sp1 3216 INTEGER :: iflavor 3217 INTEGER :: i 3218 INTEGER :: flavors 3219 CHARACTER(LEN=2) :: a 3220 CHARACTER(LEN=50) :: string 3221 DOUBLE PRECISION :: dt 3222 3223 !IF ( this%rank .NE. MOD(5,this%size)) RETURN 3224 IF ( this%rank .NE. MOD(this%size+5,this%size)) RETURN 3225 IF ( this%opt_analysis .NE. 1 ) RETURN 3226 3227 oFile = 44 3228 IF ( PRESENT(oFileIn) ) THEN 3229 oFile = oFileIn 3230 ELSE 3231 OPEN(UNIT=oFile, FILE="Correlation.dat") 3232 END IF 3233 3234 sp1 = this%samples 3235 dt = this%beta / sp1 3236 sp1 = sp1 + 1 3237 flavors = this%flavors 3238 3239 i = 3*flavors + 1 3240 WRITE(a,'(I2)') i 3241 WRITE(oFile,*) "# time (/ (segement, antiseg, correl), i=1, flavor/)" 3242 string = '(1x,'//TRIM(ADJUSTL(a))//'F19.15)' 3243 DO itime = 1, sp1 3244 WRITE(oFile,string) DBLE(itime-1)*dt, & 3245 (/ ( & 3246 (/ ( this%measCorrelation(itime, i, iflavor), i=1,3) /) & 3247 , iflavor=1, flavors) /) 3248 END DO 3249 3250 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3251 3252 END SUBROUTINE Ctqmc_printCorrelation
ABINIT/m_Ctqmc/Ctqmc_printD [ Functions ]
NAME
Ctqmc_printD
FUNCTION
print individual double occupancy
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3031 SUBROUTINE Ctqmc_printD(this,oFileIn) 3032 3033 !Arguments ------------------------------------ 3034 TYPE(Ctqmc) , INTENT(IN) :: this 3035 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3036 !Local variables ------------------------------ 3037 INTEGER :: oFile 3038 INTEGER :: iflavor1 3039 INTEGER :: iflavor2 3040 3041 !IF ( this%rank .NE. MOD(2,this%size)) RETURN 3042 IF ( this%rank .NE. MOD(this%size+2,this%size)) RETURN 3043 3044 oFile = 41 3045 IF ( PRESENT(oFileIn) ) THEN 3046 oFile = oFileIn 3047 ELSE 3048 OPEN(UNIT=oFile, FILE="D.dat") 3049 END IF 3050 3051 DO iflavor1 = 1, this%flavors 3052 DO iflavor2 = iflavor1+1, this%flavors 3053 WRITE(oFile,'(1x,A8,I4,A1,I4,A3,ES21.14)') "Orbitals", iflavor1, "-", iflavor2, " : ", this%measDE(iflavor2,iflavor1) 3054 END DO 3055 END DO 3056 3057 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3058 3059 END SUBROUTINE Ctqmc_printD
ABINIT/m_Ctqmc/Ctqmc_printE [ Functions ]
NAME
Ctqmc_printE
FUNCTION
print energy and noise
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3087 SUBROUTINE Ctqmc_printE(this,oFileIn) 3088 3089 !Arguments ------------------------------------ 3090 TYPE(Ctqmc) , INTENT(IN) :: this 3091 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3092 !Local variables ------------------------------ 3093 INTEGER :: oFile 3094 DOUBLE PRECISION :: E 3095 DOUBLE PRECISION :: Noise 3096 3097 !IF ( this%rank .NE. MOD(3,this%size)) RETURN 3098 IF ( this%rank .NE. MOD(this%size+3,this%size)) RETURN 3099 3100 oFile = 42 3101 IF ( PRESENT(oFileIn) ) THEN 3102 oFile = oFileIn 3103 ELSE 3104 OPEN(UNIT=oFile, FILE="BetaENoise.dat") 3105 END IF 3106 3107 CALL Ctqmc_getE(this,E,Noise) 3108 3109 WRITE(oFile,'(1x,F5.2,A2,ES21.14,A2,ES21.14)') this%beta, " ", E, " ", Noise 3110 3111 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3112 3113 END SUBROUTINE Ctqmc_printE
ABINIT/m_Ctqmc/Ctqmc_printGreen [ Functions ]
NAME
Ctqmc_printGreen
FUNCTION
print green functions
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2952 SUBROUTINE Ctqmc_printGreen(this, oFileIn) 2953 2954 !Arguments ------------------------------------ 2955 TYPE(Ctqmc) , INTENT(IN) :: this 2956 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 2957 !Local variables ------------------------------ 2958 INTEGER :: oFile 2959 INTEGER :: itime 2960 INTEGER :: sp1 2961 INTEGER :: iflavor 2962 INTEGER :: flavors 2963 CHARACTER(LEN=4) :: cflavors 2964 CHARACTER(LEN=50) :: string 2965 DOUBLE PRECISION :: dt 2966 DOUBLE PRECISION :: sweeps 2967 2968 !IF ( this%rank .NE. MOD(1,this%size)) RETURN 2969 IF ( this%rank .NE. MOD(this%size+1,this%size)) RETURN 2970 2971 oFile = 40 2972 IF ( PRESENT(oFileIn) ) THEN 2973 oFile = oFileIn 2974 ELSE 2975 OPEN(UNIT=oFile, FILE="Gtau.dat") 2976 END IF 2977 2978 sp1 = this%samples 2979 dt = this%beta / DBLE(sp1) 2980 sp1 = sp1 + 1 2981 flavors = this%flavors 2982 sweeps = DBLE(this%sweeps)*DBLE(this%size) 2983 2984 IF ( this%opt_noise .EQ. 1) THEN 2985 WRITE(cflavors,'(I4)') 2*flavors+1 2986 string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)' 2987 DO itime = 1, sp1 2988 WRITE(oFile,string) DBLE(itime-1)*dt, & 2989 (/ (this%Greens(iflavor)%oper(itime), iflavor=1, flavors) /), & 2990 (/ (this%abNoiseG(1,itime,iflavor)*(sweeps)**this%abNoiseG(2,itime,iflavor), iflavor=1, flavors) /) 2991 END DO 2992 ELSE 2993 WRITE(cflavors,'(I4)') flavors+1 2994 string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)' 2995 DO itime = 1, sp1 2996 WRITE(oFile,string) DBLE(itime-1)*dt, & 2997 (/ (this%Greens(iflavor)%oper(itime), iflavor=1, flavors) /) 2998 END DO 2999 END IF 3000 3001 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3002 3003 END SUBROUTINE Ctqmc_printGreen
ABINIT/m_Ctqmc/Ctqmc_printPerturbation [ Functions ]
NAME
Ctqmc_printPerturbation
FUNCTION
print perturbation order
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc oFileIn=file stream
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
3143 SUBROUTINE Ctqmc_printPerturbation(this, oFileIn) 3144 3145 !Arguments ------------------------------------ 3146 TYPE(Ctqmc) , INTENT(IN) :: this 3147 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3148 !Local variables------------------------------- 3149 INTEGER :: oFile 3150 INTEGER :: iorder 3151 INTEGER :: order 3152 INTEGER :: iflavor 3153 INTEGER :: flavors 3154 CHARACTER(LEN=2) :: a 3155 CHARACTER(LEN=50) :: string 3156 3157 !IF ( this%rank .NE. MOD(4,this%size)) RETURN 3158 IF ( this%rank .NE. MOD(this%size+4,this%size)) RETURN 3159 IF ( this%opt_order .LE. 0 ) RETURN 3160 3161 oFile = 43 3162 IF ( PRESENT(oFileIn) ) THEN 3163 oFile = oFileIn 3164 ELSE 3165 OPEN(UNIT=oFile, FILE="Perturbation.dat") 3166 END IF 3167 3168 order = this%opt_order 3169 flavors = this%flavors 3170 3171 WRITE(a,'(I2)') flavors 3172 string = '(I5,'//TRIM(ADJUSTL(a))//'F19.15)' 3173 DO iorder = 1, order 3174 WRITE(oFile,string) iorder-1, & 3175 (/ (this%measPerturbation(iorder, iflavor), iflavor=1, flavors) /) 3176 END DO 3177 3178 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3179 END SUBROUTINE Ctqmc_printPerturbation
ABINIT/m_Ctqmc/Ctqmc_printQMC [ Functions ]
NAME
Ctqmc_printQMC
FUNCTION
print ctqmc statistics
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2850 SUBROUTINE Ctqmc_printQMC(this) 2851 2852 !Arguments ------------------------------------ 2853 TYPE(Ctqmc), INTENT(INOUT) :: this 2854 !Local variables ------------------------------ 2855 INTEGER :: ostream 2856 INTEGER :: iflavor 2857 DOUBLE PRECISION :: sweeps 2858 DOUBLE PRECISION :: invSweeps 2859 CHARACTER(LEN=2) :: a 2860 CHARACTER(LEN=15) :: string 2861 2862 !IF ( this%rank .NE. 0) RETURN 2863 IF ( this%rank .NE. MOD(this%size,this%size)) RETURN 2864 2865 ostream = this%ostream 2866 sweeps = DBLE(this%sweeps) 2867 invSweeps = 1.d0/sweeps 2868 2869 WRITE(ostream,'(1x,F13.0,A11,F10.2,A12,I5,A5)') sweeps*DBLE(this%size), " sweeps in ", this%runTime, & 2870 " seconds on ", this%size, " CPUs" 2871 WRITE(ostream,'(A28,F6.2)') "Segments added [%] : ", this%stats(CTQMC_SEGME+CTQMC_ADDED)*invSweeps*100.d0 2872 WRITE(ostream,'(A28,F6.2)') "Segments removed [%] : ", this%stats(CTQMC_SEGME+CTQMC_REMOV)*invSweeps*100.d0 2873 WRITE(ostream,'(A28,F6.2)') "Segments sign [%] : ", this%stats(CTQMC_SEGME+CTQMC_DETSI)*invSweeps*100.d0 2874 WRITE(ostream,'(A28,F6.2)') "Anti-segments added [%] : ", this%stats(CTQMC_ANTIS+CTQMC_ADDED)*invSweeps*100.d0 2875 WRITE(ostream,'(A28,F6.2)') "Anti-segments removed [%] : ", this%stats(CTQMC_ANTIS+CTQMC_REMOV)*invSweeps*100.d0 2876 WRITE(ostream,'(A28,F6.2)') "Anti-segments sign [%] : ", this%stats(CTQMC_ANTIS+CTQMC_DETSI)*invSweeps*100.d0 2877 IF ( this%modGlobalMove(1) .LT. this%sweeps + 1 ) THEN 2878 WRITE(ostream,'(A28,F6.2)') "Global Move [%] : ", this%swap *invSweeps*100.d0*this%modGlobalMove(1) 2879 WRITE(ostream,'(A28,F6.2)') "Global Move Reduced [%] : ", this%swap / DBLE(this%modGlobalMove(2))*100.d0 2880 END IF 2881 !#ifdef CTCtqmc_CHECK 2882 IF ( this%opt_check .EQ. 1 .OR. this%opt_check .EQ. 3 ) & 2883 WRITE(ostream,'(A28,E22.14)') "Impurity test [%] : ", this%errorImpurity*100.d0 2884 IF ( this%opt_check .GE. 2 ) & 2885 WRITE(ostream,'(A28,E22.14)') "Bath test [%] : ", this%errorBath *100.d0 2886 !#endif 2887 WRITE(ostream,'(A28,ES22.14,A5,ES21.14)') "<Epot> [U] : ", this%measDE(1,1), " +/- ",& 2888 !#ifdef HAVE_MPI 2889 SUM(this%Impurity%mat_U)/(this%flavors*(this%flavors-1)) * this%a_Noise*(sweeps*DBLE(this%size))**this%b_Noise 2890 !#else 2891 ! this%a_Noise*(sweeps)**this%b_Noise 2892 !#endif 2893 WRITE(ostream,'(A28,F8.4,A3,F7.4)') "Noise [/U] : ", this%a_Noise, " x^", this%b_Noise 2894 WRITE(ostream,'(A28,E10.2)') "Niquist puls. [/beta] : ", ACOS(-1.d0)*this%inv_dt 2895 WRITE(ostream,'(A28,E22.14)') "Max Acc. Epot Error [U] : ", this%measDE(2,2)/(this%beta*this%modNoise1*2.d0)*sweeps 2896 2897 !WRITE(ostream,'(A28,F7.4,A3,F7.4,A4,E20.14)') "Noise [G(tau)] : ", this%a_Noise(2), "x^", this%b_Noise(2), " -> ", & 2898 !this%a_Noise(2)*(sweeps*DBLE(this%size))**this%b_Noise(2) 2899 IF ( this%opt_order .GT. 0 ) THEN 2900 WRITE(a,'(I2)') this%flavors 2901 string = '(A28,'//TRIM(ADJUSTL(a))//'(1x,I3))' 2902 WRITE(ostream,string) "Perturbation orders : ", & 2903 (/ (MAXLOC(this%measPerturbation(:, iflavor))-1, iflavor=1, this%flavors) /) 2904 END IF 2905 !CALL FLUSH(this%ostream) 2906 IF ( ABS(((this%stats(CTQMC_SEGME+CTQMC_ADDED) *invSweeps*100.d0) / & 2907 (this%stats(CTQMC_SEGME+CTQMC_REMOV) *invSweeps*100.d0) - 1.d0)) .GE. 0.02d0 & 2908 .OR. ABS(((this%stats(CTQMC_ANTIS+CTQMC_ADDED)*invSweeps*100.d0) / & 2909 (this%stats(CTQMC_ANTIS+CTQMC_REMOV)*invSweeps*100.d0) - 1.d0)) .GE. 0.02d0 ) & 2910 THEN 2911 CALL WARNALL("Ctqmc_printQMC : bad statistic according to moves. Increase sweeps") 2912 END IF 2913 ! Check sign problem for diagonal hybridization. 2914 IF ( (this%stats(CTQMC_SEGME+CTQMC_DETSI) + this%stats(CTQMC_ANTIS+CTQMC_DETSI)) .GT. 1.d-10 ) THEN 2915 CALL WARNALL("Ctqmc_printQMC : at least one negative sign occured. There might be a bug in the CT-QMC") 2916 END IF 2917 2918 IF ( ABS(this%b_Noise+0.5)/0.5d0 .GE. 0.05d0 ) & 2919 CALL WARNALL("Ctqmc_printQMC : bad statistic according to Noise. Increase sweeps") 2920 ! IF ( ISNAN(this%a_Noise) .OR. ISNAN(this%a_Noise) ) & 2921 ! CALL WARNALL("Ctqmc_printQMC : NaN appeared. Increase sweeps ") 2922 2923 2924 END SUBROUTINE Ctqmc_printQMC
ABINIT/m_Ctqmc/Ctqmc_printSpectra [ Functions ]
NAME
Ctqmc_printSpectra
FUNCTION
print fourier transform of time evolution of number of electrons
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3281 SUBROUTINE Ctqmc_printSpectra(this, oFileIn) 3282 3283 !Arguments ------------------------------------ 3284 TYPE(Ctqmc) , INTENT(IN) :: this 3285 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3286 !Local variables ------------------------------ 3287 INTEGER :: oFile 3288 INTEGER :: flavors 3289 INTEGER :: indDensity 3290 INTEGER :: endDensity 3291 CHARACTER(LEN=4) :: a 3292 CHARACTER(LEN=16) :: formatSpectra 3293 3294 !IF ( this%rank .NE. MOD(6,this%size)) RETURN 3295 IF ( this%opt_spectra .LT. 1 ) RETURN 3296 3297 oFile = 45+this%rank 3298 a ="0000" 3299 WRITE(a,'(I4)') this%rank 3300 IF ( PRESENT(oFileIn) ) THEN 3301 oFile = oFileIn 3302 ELSE 3303 OPEN(UNIT=oFile, FILE="Markov_"//TRIM(ADJUSTL(a))//".dat") 3304 END IF 3305 3306 flavors = this%flavors 3307 WRITE(a,'(I4)') flavors+1 3308 formatSpectra ='(1x,'//TRIM(ADJUSTL(a))//'ES22.14)' 3309 WRITE(oFile,*) "# freq[/hermalization] FFT" 3310 3311 endDensity = SIZE(this%density,2) 3312 DO WHILE ( this%density(flavors+1,endDensity) .EQ. -1 ) 3313 endDensity = endDensity -1 3314 END DO 3315 3316 DO indDensity = 1, endDensity 3317 WRITE(oFile,formatSpectra) this%density(flavors+1,indDensity), this%density(1:flavors,indDensity) 3318 END DO 3319 3320 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3321 3322 END SUBROUTINE Ctqmc_printSpectra
ABINIT/m_Ctqmc/Ctqmc_reset [ Functions ]
NAME
Ctqmc_reset
FUNCTION
reset a ctqmc simulation
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1012 SUBROUTINE Ctqmc_reset(this) 1013 1014 !Arguments ------------------------------------ 1015 TYPE(Ctqmc), INTENT(INOUT) :: this 1016 !Local variables ------------------------------ 1017 INTEGER :: iflavor 1018 DOUBLE PRECISION :: sweeps 1019 1020 DO iflavor = 1, this%flavors 1021 CALL GreenHyb_reset(this%Greens(iflavor)) 1022 END DO 1023 CALL Ctqmc_clear(this) 1024 CALL ImpurityOperator_reset(this%Impurity) 1025 CALL BathOperator_reset (this%Bath) 1026 this%measN(3,:) = 0.d0 1027 !complete restart -> measN=0 1028 this%done = .FALSE. 1029 this%set = .FALSE. 1030 this%inF = .FALSE. 1031 this%opt_movie = 0 1032 this%opt_analysis = 0 1033 this%opt_order = 0 1034 this%opt_check = 0 1035 this%opt_noise = 0 1036 this%opt_spectra = 0 1037 this%opt_levels = 0 1038 sweeps = DBLE(this%sweeps)*DBLE(this%size) 1039 CALL Ctqmc_setSweeps(this, sweeps) 1040 !#ifdef HAVE_MPI 1041 ! CALL MPI_BARRIER(this%MY_COMM,iflavor) 1042 ! IF ( this%rank .EQ. 0 ) & 1043 !#endif 1044 ! WRITE(this%ostream,'(A9)') "QMC reset" 1045 ! CALL FLUSH(this%ostream) 1046 END SUBROUTINE Ctqmc_reset
ABINIT/m_Ctqmc/Ctqmc_run [ Functions ]
NAME
Ctqmc_run
FUNCTION
set all options and run a simulation
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc opt_order=maximal perturbation order to scope opt_movie=draw a movie of the simulation opt_analysis=compute correlation functions opt_check=check fast calculations opt_noise=compute noise for green function opt_spectra=fourier transform of the time evolution of the number of electrons opt_gMove=steps without global move
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1305 SUBROUTINE Ctqmc_run(this,opt_order,opt_histo,opt_movie,opt_analysis,opt_check,opt_noise,opt_spectra,opt_gMove) 1306 1307 1308 #ifdef HAVE_MPI1 1309 include 'mpif.h' 1310 #endif 1311 !Arguments ------------------------------------ 1312 TYPE(Ctqmc), INTENT(INOUT) :: this 1313 INTEGER, OPTIONAL, INTENT(IN ) :: opt_order 1314 INTEGER, OPTIONAL, INTENT(IN ) :: opt_histo 1315 INTEGER, OPTIONAL, INTENT(IN ) :: opt_movie 1316 INTEGER, OPTIONAL, INTENT(IN ) :: opt_analysis 1317 INTEGER, OPTIONAL, INTENT(IN ) :: opt_check 1318 INTEGER, OPTIONAL, INTENT(IN ) :: opt_noise 1319 INTEGER, OPTIONAL, INTENT(IN ) :: opt_spectra 1320 INTEGER, OPTIONAL, INTENT(IN ) :: opt_gMove 1321 !Local variables ------------------------------ 1322 #ifdef HAVE_MPI 1323 INTEGER :: ierr 1324 #endif 1325 !#ifdef CTCtqmc_MOVIE 1326 INTEGER :: ilatex 1327 CHARACTER(LEN=4) :: Cchar 1328 !#endif 1329 DOUBLE PRECISION :: estimatedTime 1330 1331 IF ( .NOT. this%set ) & 1332 CALL ERROR("Ctqmc_run : QMC not set up ") 1333 IF ( .NOT. this%setU ) & 1334 CALL ERROR("Ctqmc_run : QMC does not have a U this ") 1335 1336 1337 ! OPTIONS of the run 1338 IF ( PRESENT( opt_check ) ) THEN 1339 this%opt_check = opt_check 1340 CALL ImpurityOperator_doCheck(this%Impurity,opt_check) 1341 CALL BathOperator_doCheck(this%Bath,opt_check) 1342 END IF 1343 IF ( PRESENT( opt_movie ) ) & 1344 this%opt_movie = opt_movie 1345 IF ( PRESENT( opt_analysis ) ) & 1346 this%opt_analysis = opt_analysis 1347 IF ( PRESENT ( opt_order ) ) & 1348 this%opt_order = opt_order 1349 IF ( PRESENT ( opt_histo ) ) & 1350 this%opt_histo = opt_histo 1351 IF ( PRESENT ( opt_noise ) ) THEN 1352 this%opt_noise = opt_noise 1353 END IF 1354 IF ( PRESENT ( opt_spectra ) ) & 1355 this%opt_spectra = opt_spectra 1356 1357 this%modGlobalMove(1) = this%sweeps+1 ! No Global Move 1358 this%modGlobalMove(2) = 0 1359 IF ( PRESENT ( opt_gMove ) ) THEN 1360 IF ( opt_gMove .LE. 0 .OR. opt_gMove .GT. this%sweeps ) THEN 1361 this%modGlobalMove(1) = this%sweeps+1 1362 CALL WARNALL("Ctqmc_run : global moves option is <= 0 or > sweeps/cpu -> No global Moves") 1363 ELSE 1364 this%modGlobalMove(1) = opt_gMove 1365 END IF 1366 END IF 1367 1368 CALL Ctqmc_allocateOpt(this) 1369 1370 !#ifdef CTCtqmc_MOVIE 1371 ilatex = 0 1372 IF ( this%opt_movie .EQ. 1 ) THEN 1373 Cchar ="0000" 1374 WRITE(Cchar,'(I4)') this%rank 1375 ilatex = 87+this%rank 1376 OPEN(UNIT=ilatex, FILE="Movie_"//TRIM(ADJUSTL(Cchar))//".tex") 1377 WRITE(ilatex,'(A)') "\documentclass{beamer}" 1378 WRITE(ilatex,'(A)') "\usepackage{color}" 1379 WRITE(ilatex,'(A)') "\setbeamersize{sidebar width left=0pt}" 1380 WRITE(ilatex,'(A)') "\setbeamersize{sidebar width right=0pt}" 1381 WRITE(ilatex,'(A)') "\setbeamersize{text width left=0pt}" 1382 WRITE(ilatex,'(A)') "\setbeamersize{text width right=0pt}" 1383 WRITE(ilatex,*) 1384 WRITE(ilatex,'(A)') "\begin{document}" 1385 WRITE(ilatex,*) 1386 END IF 1387 !#endif 1388 1389 IF ( this%rank .EQ. 0 ) THEN 1390 WRITE(this%ostream,'(A29)') "Starting QMC (Thermalization)" 1391 END IF 1392 1393 !================================= 1394 ! STARTING THERMALIZATION 1395 !================================= 1396 CALL Ctqmc_loop(this,this%thermalization,ilatex) 1397 !================================= 1398 ! ENDING THERMALIZATION 1399 !================================= 1400 1401 estimatedTime = this%runTime 1402 #ifdef HAVE_MPI 1403 CALL MPI_REDUCE(this%runTime, estimatedTime, 1, MPI_DOUBLE_PRECISION, MPI_MAX, & 1404 0, this%MY_COMM, ierr) 1405 #endif 1406 1407 IF ( this%rank .EQ. 0 ) THEN 1408 WRITE(this%ostream,'(A26,I6,A11)') "Thermalization done in ", CEILING(estimatedTime), " seconds" 1409 WRITE(this%ostream,'(A25,I7,A15,I5,A5)') "The QMC should run in ", & 1410 CEILING(estimatedTime*DBLE(this%sweeps)/DBLE(this%thermalization)),& 1411 " seconds on ", this%size, " CPUs" 1412 END IF 1413 1414 !================================= 1415 ! CLEANING CTQMC 1416 !================================= 1417 CALL Ctqmc_clear(this) 1418 1419 !================================= 1420 ! STARTING CTQMC 1421 !================================= 1422 CALL Ctqmc_loop(this,this%sweeps,ilatex) 1423 !================================= 1424 ! ENDING CTQMC 1425 !================================= 1426 1427 IF ( this%opt_movie .EQ. 1 ) THEN 1428 WRITE(ilatex,*) "" 1429 WRITE(ilatex,'(A14)') "\end{document}" 1430 CLOSE(ilatex) 1431 END IF 1432 1433 this%done = .TRUE. 1434 1435 END SUBROUTINE Ctqmc_run
ABINIT/m_Ctqmc/Ctqmc_setG0wTab [ Functions ]
NAME
Ctqmc_setG0wTab
FUNCTION
Set Gow from input array
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc Gomega=G0w opt_fk=F is already inversed with out iwn
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
814 SUBROUTINE Ctqmc_setG0wTab(this,Gomega,opt_fk) 815 816 !Arguments ------------------------------------ 817 TYPE(Ctqmc), INTENT(INOUT) :: this 818 COMPLEX(KIND=8), DIMENSION(:,:), INTENT(IN ) :: Gomega 819 INTEGER , INTENT(IN ) :: opt_fk 820 !Local variable ------------------------------- 821 DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: F 822 823 IF ( .NOT. this%para ) & 824 CALL ERROR("Ctqmc_setG0wTab : Ctqmc_setParameters never called ") 825 826 MALLOC(F,(1:this%samples+1,1:this%flavors)) 827 CALL Ctqmc_computeF(this,Gomega, F, opt_fk) ! mu is changed 828 CALL BathOperator_setF(this%Bath, F) 829 !CALL BathOperator_printF(this%Bath) 830 FREE(F) 831 832 IF ( this%opt_levels .NE. 1 ) THEN ! We compute the mu by hand in computeF 833 CALL ImpurityOperator_setMu(this%Impurity,this%mu) 834 END IF 835 836 this%inF = .TRUE. 837 this%set = .TRUE. 838 839 END SUBROUTINE Ctqmc_setG0wTab
ABINIT/m_Ctqmc/Ctqmc_setMu [ Functions ]
NAME
Ctqmc_setMu
FUNCTION
impose energy levels
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc levels=energy levels vector
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
1075 SUBROUTINE Ctqmc_setMu(this, levels) 1076 1077 !Arguments ------------------------------------ 1078 TYPE(Ctqmc) , INTENT(INOUT) :: this 1079 DOUBLE PRECISION, DIMENSION(:), INTENT(IN ) :: levels 1080 1081 IF ( this%flavors .NE. SIZE(levels,1) ) & 1082 CALL WARNALL("Ctqmc_setMu : Taking energy levels from weiss G(iw)") 1083 1084 this%mu(:)=-levels ! levels = \epsilon_j - \mu 1085 !this%mu =\tilde{\mu} = \mu -\epsilon_j 1086 CALL ImpurityOperator_setMu(this%Impurity,this%mu) 1087 this%opt_levels = 1 1088 END SUBROUTINE Ctqmc_setMu
ABINIT/m_Ctqmc/Ctqmc_setParameters [ Functions ]
NAME
Ctqmc_setParameters
FUNCTION
set all parameters and operators
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc buffer=input parameters
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
456 SUBROUTINE Ctqmc_setParameters(this,buffer) 457 458 !Arguments ------------------------------------ 459 TYPE(Ctqmc), INTENT(INOUT) :: this 460 DOUBLE PRECISION, DIMENSION(1:9), INTENT(IN ) :: buffer 461 462 463 this%thermalization = INT(buffer(3)) !this%thermalization 464 CALL Ctqmc_setSeed(this,INT(buffer(1))) 465 CALL Ctqmc_setSweeps(this,buffer(2)) 466 467 this%measurements = INT(buffer(4)) !this%measurements 468 this%flavors = INT(buffer(5)) !this%flavors 469 this%samples = INT(buffer(6)) !this%samples 470 this%beta = buffer(7) !this%beta 471 this%U = buffer(8) !U 472 ! this%mu = buffer(9) !this%mu 473 !this%Wmax = INT(buffer(9)) !Freq 474 !#ifdef CTCtqmc_ANALYSIS 475 ! this%order = INT(buffer(10)) ! order 476 this%inv_dt = this%samples / this%beta 477 !#endif 478 479 !CALL ImpurityOperator_init(this%Impurity,this%flavors,this%beta, this%samples) 480 CALL ImpurityOperator_init(this%Impurity,this%flavors,this%beta) 481 IF ( this%U .GE. 0.d0 ) THEN 482 CALL ImpurityOperator_computeU(this%Impurity,this%U,0.d0) 483 this%setU = .TRUE. 484 END IF 485 ! this%mu = this%mu + this%Impurity%shift_mu 486 487 CALL BathOperator_init(this%Bath, this%flavors, this%samples, this%beta, INT(buffer(9))) 488 489 this%para = .TRUE. 490 491 END SUBROUTINE Ctqmc_setParameters
ABINIT/m_Ctqmc/Ctqmc_setSeed [ Functions ]
NAME
Ctqmc_setSeed
FUNCTION
initialize random number generator
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc iseed=seed from imput
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
574 SUBROUTINE Ctqmc_setSeed(this,iseed) 575 576 !Arguments ------------------------------------ 577 TYPE(Ctqmc), INTENT(INOUT) :: this 578 INTEGER , INTENT(IN ) :: iseed 579 !Local variables ------------------------------ 580 !INTEGER :: n 581 !INTEGER :: i 582 !INTEGER, DIMENSION(:), ALLOCATABLE :: seed 583 584 585 !CALL RANDOM_SEED(size = n) 586 !MALLOC(seed,(n)) 587 !seed = iseed + (/ (i - 1, i = 1, n) /) 588 589 !CALL RANDOM_SEED(PUT = seed+this%rank) 590 591 !FREE(seed) 592 593 this%seed=INT(iseed+this%rank,8) 594 595 END SUBROUTINE Ctqmc_setSeed
ABINIT/m_Ctqmc/Ctqmc_setSweeps [ Functions ]
NAME
Ctqmc_setSweeps
FUNCTION
set the number of sweeps
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc sweeps=asked sweeps
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
519 SUBROUTINE Ctqmc_setSweeps(this,sweeps) 520 521 !Arguments ------------------------------------ 522 TYPE(Ctqmc) , INTENT(INOUT) :: this 523 DOUBLE PRECISION , INTENT(IN ) :: sweeps 524 525 this%sweeps = NINT(sweeps / DBLE(this%size)) 526 ! write(6,*) this%sweeps,NINT(sweeps / DBLE(this%size)),ANINT(sweeps/DBLE(this%size)) 527 IF ( DBLE(this%sweeps) .NE. ANINT(sweeps/DBLE(this%size)) ) & 528 CALL ERROR("Ctqmc_setSweeps : sweeps is negative or too big ") 529 IF ( this%sweeps .LT. 2*CTQMC_SLICE1 ) THEN !202 530 CALL WARNALL("Ctqmc_setSweeps : # sweeps automtically changed ") 531 this%sweeps = 2*CTQMC_SLICE1 532 ! ELSE IF ( this%sweeps .LT. this%thermalization ) THEN 533 ! CALL WARNALL("Ctqmc_setSweeps : Thermalization > sweeps / cpu -> auto fix") 534 ! this%sweeps = this%thermalization 535 END IF 536 IF ( DBLE(NINT(DBLE(this%sweeps)*DBLE(this%size)/DBLE(CTQMC_SLICE1))) .NE. & 537 ANINT(DBLE(this%sweeps)*DBLE(this%size)/DBLE(CTQMC_SLICE1)) ) THEN 538 this%modNoise1 = this%sweeps 539 ELSE 540 this%modNoise1 = MIN(this%sweeps,INT(DBLE(this%sweeps)*DBLE(this%size) / DBLE(CTQMC_SLICE1))) !101 541 END IF 542 this%modNoise2 = MAX(this%modNoise1 / CTQMC_SLICE2, 1) ! 100 543 ! this%modGlobalMove(1) = this%thermalization / 10 + 1 544 ! this%modGlobalMove(2) = 0 545 546 END SUBROUTINE Ctqmc_setSweeps
ABINIT/m_Ctqmc/Ctqmc_setU [ Functions ]
NAME
Ctqmc_setU
FUNCTION
set the interaction this
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc matU=interaction this
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
898 SUBROUTINE Ctqmc_setU(this,matU) 899 900 !Arguments ------------------------------------ 901 TYPE(Ctqmc), INTENT(INOUT) :: this 902 !Local variables ------------------------------ 903 DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: matU 904 905 IF ( SIZE(matU) .NE. this%flavors*this%flavors ) & 906 CALL ERROR("Ctqmc_setU : Wrong interaction this (size) ") 907 908 CALL ImpurityOperator_setUmat(this%Impurity, matU) 909 this%setU = .TRUE. 910 END SUBROUTINE Ctqmc_setU
ABINIT/m_Ctqmc/Ctqmc_symmetrizeGreen [ Functions ]
NAME
Ctqmc_symmetrizeGreen
FUNCTION
optionnaly symmetrize the green functions
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc syms=weight factors
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
2564 SUBROUTINE Ctqmc_symmetrizeGreen(this, syms) 2565 2566 !Arguments ------------------------------------ 2567 TYPE(Ctqmc) , INTENT(INOUT) :: this 2568 DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN ) :: syms 2569 !Local variables ------------------------------ 2570 INTEGER :: iflavor1 2571 INTEGER :: iflavor2 2572 INTEGER :: flavors 2573 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: green_tmp 2574 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(: ) :: n_tmp 2575 2576 flavors = this%flavors 2577 IF ( SIZE(syms,1) .NE. flavors .OR. SIZE(syms,2) .NE. flavors ) THEN 2578 CALL WARNALL("Ctqmc_symmetrizeGreen : wrong opt_sym -> not symmetrizing") 2579 RETURN 2580 END IF 2581 2582 MALLOC(green_tmp,(1:this%samples+1,flavors)) 2583 green_tmp(:,:) = 0.d0 2584 MALLOC(n_tmp,(1:flavors)) 2585 n_tmp(:) = 0.d0 2586 DO iflavor1=1, flavors 2587 DO iflavor2=1,flavors 2588 green_tmp(:,iflavor1) = green_tmp(:,iflavor1) & 2589 + syms(iflavor2,iflavor1) * this%Greens(iflavor2)%oper(:) 2590 n_tmp(iflavor1) = n_tmp(iflavor1) & 2591 + syms(iflavor2,iflavor1) * this%measN(1,iflavor2) 2592 END DO 2593 END DO 2594 DO iflavor1=1, flavors 2595 this%Greens(iflavor1)%oper(:) = green_tmp(:,iflavor1) 2596 this%measN(1,iflavor1) = n_tmp(iflavor1) 2597 END DO 2598 FREE(green_tmp) 2599 FREE(n_tmp) 2600 END SUBROUTINE Ctqmc_symmetrizeGreen
ABINIT/m_Ctqmc/Ctqmc_tryAddRemove [ Functions ]
NAME
Ctqmc_tryAddRemove
FUNCTION
Try to add or remove a segment and an anti-segment
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
updated=something changed
SIDE EFFECTS
NOTES
SOURCE
1700 SUBROUTINE Ctqmc_tryAddRemove(this,updated) 1701 1702 !Arguments ------------------------------------ 1703 TYPE(Ctqmc) , INTENT(INOUT) :: this 1704 ! TYPE(BathOperator) , INTENT(INOUT) :: Bath 1705 ! TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity 1706 LOGICAL , INTENT( OUT) :: updated 1707 !Local variables ------------------------------ 1708 INTEGER :: position 1709 INTEGER , DIMENSION(1:2) :: nature ! -2 for antiseg and 1 for seg 1710 INTEGER :: i! -2 for antiseg and 1 for seg 1711 DOUBLE PRECISION :: action 1712 DOUBLE PRECISION :: beta 1713 DOUBLE PRECISION :: time1 1714 DOUBLE PRECISION :: time2 1715 DOUBLE PRECISION :: time_avail 1716 DOUBLE PRECISION :: det_ratio 1717 DOUBLE PRECISION :: imp_trace 1718 DOUBLE PRECISION :: signe 1719 DOUBLE PRECISION :: tail 1720 DOUBLE PRECISION, DIMENSION(1:2) :: CdagC_1 1721 1722 IF ( .NOT. this%set ) & 1723 CALL ERROR("Ctqmc_trySegment : QMC not set ") 1724 1725 nature(1) = CTQMC_SEGME 1726 nature(2) = CTQMC_ANTIS 1727 beta = this%beta 1728 1729 updated = .FALSE. 1730 tail = DBLE(this%Impurity%particles(this%Impurity%activeFlavor)%tail) 1731 1732 1733 DO i = 1, 2 1734 signe = SIGN(1.d0,DBLE(nature(i))) 1735 1736 !CALL RANDOM_NUMBER(action) 1737 CALL OurRng(this%seed,action) 1738 1739 IF ( action .LT. .5d0 ) THEN ! Ajout de segment 1740 !CALL RANDOM_NUMBER(time1) 1741 CALL OurRng(this%seed,time1) 1742 time1 = time1 * beta 1743 time_avail = ImpurityOperator_getAvailableTime(this%Impurity,time1,position) * signe 1744 IF ( time_avail .GT. 0.d0 ) THEN 1745 !CALL RANDOM_NUMBER(time2) 1746 CALL OurRng(this%seed,time2) 1747 IF ( time2 .EQ. 0.d0 ) THEN 1748 CALL OurRng(this%seed,time2) ! Prevent null segment 1749 END IF 1750 time2 = time1 + time2 * time_avail 1751 CdagC_1(Cdag_) = ((1.d0+signe)*time1+(1.d0-signe)*time2)*0.5d0 1752 CdagC_1(C_ ) = ((1.d0+signe)*time2+(1.d0-signe)*time1)*0.5d0 1753 det_ratio = BathOperator_getDetAdd(this%Bath,CdagC_1,position,this%Impurity%particles(this%Impurity%activeFlavor)) 1754 imp_trace = ImpurityOperator_getTraceAdd(this%Impurity,CdagC_1) 1755 !CALL RANDOM_NUMBER(time1) 1756 CALL OurRng(this%seed,time1) 1757 IF ( det_ratio*imp_trace .LT. 0.d0 ) THEN 1758 this%stats(nature(i)+CTQMC_DETSI) = this%stats(nature(i)+CTQMC_DETSI) + 1.d0 1759 END IF 1760 IF ( (time1 * (tail + 1.d0 )) & 1761 .LT. (beta * time_avail * det_ratio * imp_trace ) ) THEN 1762 CALL ImpurityOperator_add(this%Impurity,CdagC_1,position) 1763 CALL BathOperator_setMAdd(this%bath,this%Impurity%particles(this%Impurity%activeFlavor)) 1764 this%stats(nature(i)+CTQMC_ADDED) = this%stats(nature(i)+CTQMC_ADDED) + 1.d0 1765 updated = .TRUE. .OR. updated 1766 tail = tail + 1.d0 1767 END IF 1768 END IF 1769 1770 ELSE ! Supprimer un segment 1771 IF ( tail .GT. 0.d0 ) THEN 1772 !CALL RANDOM_NUMBER(time1) 1773 CALL OurRng(this%seed,time1) 1774 position = INT(((time1 * tail) + 1.d0) * signe ) 1775 time_avail = ImpurityOperator_getAvailedTime(this%Impurity,position) 1776 det_ratio = BathOperator_getDetRemove(this%Bath,position) 1777 imp_trace = ImpurityOperator_getTraceRemove(this%Impurity,position) 1778 !CALL RANDOM_NUMBER(time1) 1779 CALL OurRng(this%seed,time1) 1780 IF ( det_ratio * imp_trace .LT. 0.d0 ) THEN 1781 this%stats(nature(i)+CTQMC_DETSI) = this%stats(nature(i)+CTQMC_DETSI) + 1.d0 1782 END IF 1783 IF ( (time1 * beta * time_avail ) & 1784 .LT. (tail * det_ratio * imp_trace) ) THEN 1785 CALL ImpurityOperator_remove(this%Impurity,position) 1786 CALL BathOperator_setMRemove(this%Bath,this%Impurity%particles(this%Impurity%activeFlavor)) 1787 this%stats(nature(i)+CTQMC_REMOV) = this%stats(nature(i)+CTQMC_REMOV) + 1.d0 1788 updated = .TRUE. .OR. updated 1789 tail = tail -1.d0 1790 END IF 1791 END IF 1792 END IF 1793 END DO 1794 END SUBROUTINE Ctqmc_tryAddRemove
ABINIT/m_Ctqmc/Ctqmc_trySwap [ Functions ]
NAME
Ctqmc_trySwap
FUNCTION
try a global move (swap to flavors)
COPYRIGHT
Copyright (C) 2013-2022 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=ctqmc
OUTPUT
flav_i=first flavor swaped flav_j=second flavor swaped
SIDE EFFECTS
NOTES
SOURCE
1823 SUBROUTINE Ctqmc_trySwap(this,flav_i,flav_j) 1824 1825 !Arguments ------------------------------------ 1826 TYPE(Ctqmc) , INTENT(INOUT) :: this 1827 ! TYPE(BathOperator) , INTENT(INOUT) :: Bath 1828 ! TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity 1829 INTEGER , INTENT( OUT) :: flav_i 1830 INTEGER , INTENT( OUT) :: flav_j 1831 !Local variables ------------------------------ 1832 INTEGER :: flavor_i 1833 INTEGER :: flavor_j 1834 DOUBLE PRECISION :: rnd 1835 DOUBLE PRECISION :: lengthi 1836 DOUBLE PRECISION :: lengthj 1837 DOUBLE PRECISION :: overlapic1 1838 DOUBLE PRECISION :: overlapjc1 1839 DOUBLE PRECISION :: overlapic2 1840 DOUBLE PRECISION :: overlapjc2 1841 DOUBLE PRECISION :: detic1 1842 DOUBLE PRECISION :: detjc1 1843 DOUBLE PRECISION :: detic2 1844 DOUBLE PRECISION :: detjc2 1845 DOUBLE PRECISION :: det_ratio 1846 DOUBLE PRECISION :: local_ratio 1847 1848 !CALL RANDOM_NUMBER(rnd) 1849 CALL OurRng(this%seed,rnd) 1850 flavor_i = NINT(rnd*DBLE(this%flavors-1.d0))+1 1851 !CALL RANDOM_NUMBER(rnd) 1852 CALL OurRng(this%seed,rnd) 1853 flavor_j = NINT(rnd*DBLE(this%flavors-1.d0))+1 1854 1855 flav_i = 0 1856 flav_j = 0 1857 1858 IF ( flavor_i .NE. flavor_j ) THEN 1859 ! On tente d'intervertir i et j 1860 ! Configuration actuelle : 1861 this%modGlobalMove(2) = this%modGlobalMove(2)+1 1862 detic1 = BathOperator_getDetF(this%Bath,flavor_i) 1863 detjc1 = BathOperator_getDetF(this%Bath,flavor_j) 1864 lengthi = ImpurityOperator_measN(this%Impurity,flavor_i) 1865 lengthj = ImpurityOperator_measN(this%Impurity,flavor_j) 1866 overlapic1 = ImpurityOperator_overlapFlavor(this%Impurity,flavor_i) 1867 overlapjc1 = ImpurityOperator_overlapFlavor(this%Impurity,flavor_j) 1868 ! Configuration nouvelle : 1869 detic2 = BathOperator_getDetF(this%Bath,flavor_i,this%Impurity%particles(flavor_j)) 1870 detjc2 = BathOperator_getDetF(this%Bath,flavor_j,this%Impurity%particles(flavor_i)) 1871 ! lengths unchanged 1872 overlapic2 = ImpurityOperator_overlapSwap(this%Impurity,flavor_i,flavor_j) 1873 overlapjc2 = ImpurityOperator_overlapSwap(this%Impurity,flavor_j,flavor_i) 1874 1875 ! IF ( detic1*detjc1 .EQ. detic2*detjc2 ) THEN 1876 ! det_ratio = 1.d0 1877 ! ELSE IF ( detic1*detjc1 .EQ. 0.d0 ) THEN 1878 ! det_ratio = detic2*detjc2 ! evite de diviser par 0 si pas de segment 1879 ! ELSE 1880 det_ratio = detic2*detjc2/(detic1*detjc1) 1881 ! END IF 1882 local_ratio = DEXP(-overlapic2*overlapjc2+overlapic1*overlapjc1 & 1883 +(lengthj-lengthi)*(this%mu(flavor_i)-this%mu(flavor_j))) 1884 ! Wloc = exp(muN-Uo) 1885 !CALL RANDOM_NUMBER(rnd) 1886 CALL OurRng(this%seed,rnd) 1887 IF ( rnd .LT. local_ratio*det_ratio ) THEN ! swap accepted 1888 CALL ImpurityOperator_swap(this%Impurity, flavor_i,flavor_j) 1889 CALL BathOperator_swap (this%Bath , flavor_i,flavor_j) 1890 this%swap = this%swap + 1.d0 1891 flav_i = flavor_i 1892 flav_j = flavor_j 1893 ! ELSE 1894 ! CALL WARN("Swap refused") 1895 ! WRITE(this%ostream,'(6E24.14)') local_ratio, det_ratio, detic1, detjc1, detic2, detjc2 1896 END IF 1897 END IF 1898 1899 END SUBROUTINE Ctqmc_trySwap
m_Ctqmc/Ctqmc [ Types ]
NAME
Ctqmc
FUNCTION
This structured datatype contains the necessary data
COPYRIGHT
Copyright (C) 2013-2022 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
73 TYPE, PUBLIC :: Ctqmc 74 75 LOGICAL _PRIVATE :: init = .FALSE. 76 ! Flag: is MC initialized 77 78 LOGICAL _PRIVATE :: set = .FALSE. 79 ! Flag: ?? 80 81 LOGICAL _PRIVATE :: setU = .FALSE. 82 ! Flag: is U Set ? 83 84 LOGICAL _PRIVATE :: inF = .FALSE. 85 ! Flag: is hybridization fct in input ? 86 87 LOGICAL _PRIVATE :: done = .FALSE. 88 ! Flag: is MC terminated ? 89 90 LOGICAL _PRIVATE :: para = .FALSE. 91 ! Flag: do we have parameters in input 92 93 LOGICAL _PRIVATE :: have_MPI = .FALSE. 94 ! Flag: 95 96 INTEGER _PRIVATE :: opt_movie = 0 97 ! 98 99 INTEGER _PRIVATE :: opt_analysis = 0 100 ! correlations 101 102 INTEGER _PRIVATE :: opt_check = 0 103 ! various check 0 104 ! various check 1 impurity 105 ! various check 2 bath 106 ! various check 3 both 107 108 INTEGER _PRIVATE :: opt_order = 0 109 ! nb of segments max for analysis 110 111 INTEGER _PRIVATE :: opt_histo = 0 112 ! Enable histograms 113 114 INTEGER _PRIVATE :: opt_noise = 0 115 ! compute noise 116 117 INTEGER _PRIVATE :: opt_spectra = 0 118 ! markov chain FT (correlation time) 119 120 INTEGER _PRIVATE :: opt_levels = 0 121 ! do we have energy levels 122 123 INTEGER _PRIVATE :: flavors 124 ! 125 126 INTEGER _PRIVATE :: measurements 127 ! nb of measure in the MC 128 129 INTEGER _PRIVATE :: samples 130 ! nb of L points 131 132 INTEGER(8) _PRIVATE :: seed 133 ! 134 135 INTEGER _PRIVATE :: sweeps 136 ! 137 138 INTEGER _PRIVATE :: thermalization 139 ! 140 141 INTEGER _PRIVATE :: ostream 142 ! output file 143 144 INTEGER _PRIVATE :: istream 145 ! input file 146 147 INTEGER _PRIVATE :: modNoise1 148 ! measure the noise each modNoise1 149 150 INTEGER _PRIVATE :: modNoise2 151 ! measure the noise each modNoise2 152 153 INTEGER _PRIVATE :: activeFlavor 154 ! orbital on which one do sth now 155 156 INTEGER, DIMENSION(1:2) _PRIVATE :: modGlobalMove 157 ! 1: gloabl move each modglobalmove(1) 158 ! 2: we have done modglobalmove(2) for two different orbitals. 159 160 INTEGER _PRIVATE :: Wmax 161 ! Max freq for FT 162 163 DOUBLE PRECISION, DIMENSION(1:6) _PRIVATE :: stats 164 ! to now how many negative determinant, antisegments,seeme.e.twfs...j 165 166 DOUBLE PRECISION _PRIVATE :: swap 167 ! nb of successfull GM 168 169 INTEGER _PRIVATE :: MY_COMM 170 ! 171 172 INTEGER _PRIVATE :: rank 173 ! 174 175 INTEGER _PRIVATE :: size 176 ! size of MY_COMM 177 178 DOUBLE PRECISION _PRIVATE :: runTime ! time for the run routine 179 ! 180 181 DOUBLE PRECISION _PRIVATE :: beta 182 ! 183 184 DOUBLE PRECISION _PRIVATE :: U 185 186 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) _PRIVATE :: mu 187 ! levels 188 189 TYPE(GreenHyb) , ALLOCATABLE, DIMENSION(: ) _PRIVATE :: Greens 190 191 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) _PRIVATE :: measN 192 ! measure of occupations (3or4,flavor) 193 194 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) _PRIVATE :: measDE 195 ! (flavor,flavor) double occupancies 196 ! (1,1): total energy of correlation. 197 198 DOUBLE PRECISION _PRIVATE :: a_Noise 199 ! Noise a exp (-bx) for the noise 200 201 DOUBLE PRECISION _PRIVATE :: b_Noise 202 ! Noise a exp (-bx) for the noise 203 204 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: abNoiseG !(ab,tau,flavor) 205 ! Noise but for G 206 207 TYPE(Vector) , DIMENSION(1:2) _PRIVATE :: measNoise 208 TYPE(Vector), ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: measNoiseG !(tau,flavor,mod) 209 ! accumulate each value relataed to measurenoise 1 2 210 211 DOUBLE PRECISION _PRIVATE :: inv_dt 212 ! 1/(beta/L) 213 214 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) _PRIVATE :: measPerturbation 215 ! opt_order,nflavor 216 217 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: occup_histo_time 218 ! nflavor 219 220 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: measCorrelation 221 ! segment,antisegment,nflavor,nflavor 222 223 DOUBLE PRECISION _PRIVATE :: errorImpurity 224 ! check 225 226 DOUBLE PRECISION _PRIVATE :: errorBath 227 ! for check 228 229 TYPE(BathOperator) _PRIVATE :: Bath 230 231 TYPE(ImpurityOperator) _PRIVATE :: Impurity 232 233 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) _PRIVATE :: density 234 235 END TYPE Ctqmc