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
- ABINIT/m_Ctqmcoffdiag/Ctqmc_setMagmom
- 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
NOTES
SOURCE
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 use m_io_tools, only : open_file 38 use defs_basis 39 #ifdef HAVE_MPI2 40 USE mpi 41 #endif 42 43 IMPLICIT NONE
ABINIT/m_Ctqmc/Ctqmc_allocateAll [ Functions ]
NAME
Ctqmc_allocateAll
FUNCTION
Allocate all non option varibales
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
642 SUBROUTINE Ctqmc_allocateAll(this) 643 644 !Arguments ------------------------------------ 645 TYPE(Ctqmc), INTENT(INOUT) :: this 646 !Local variables ------------------------------ 647 INTEGER :: flavors 648 649 IF ( .NOT. this%para ) & 650 CALL ERROR("Ctqmc_allocateAll : Ctqmc_setParameters never called ") 651 652 flavors = this%flavors 653 654 DT_FREEIF(this%Greens) 655 DT_MALLOC(this%Greens,(1:flavors)) 656 657 FREEIF(this%measN) 658 MALLOC(this%measN,(1:4,1:flavors)) 659 this%measN = 0.d0 660 661 FREEIF(this%measDE) 662 MALLOC(this%measDE,(1:flavors,1:flavors) ) 663 this%measDE = 0.d0 664 665 FREEIF(this%mu) 666 #ifdef FC_LLVM 667 ! LLVM 16 doesn't recognize this macro here 668 MALLOC(this%mu, (1:flavors) ) 669 #else 670 MALLOC(this%mu, (1:flavors)) 671 #endif 672 this%mu = 0.d0 673 END SUBROUTINE Ctqmc_allocateAll
ABINIT/m_Ctqmc/Ctqmc_allocateOpt [ Functions ]
NAME
Ctqmc_allocateOpt
FUNCTION
allocate all option variables
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
701 SUBROUTINE Ctqmc_allocateOpt(this) 702 703 !Arguments ------------------------------------ 704 TYPE(Ctqmc), INTENT(INOUT) :: this 705 !Local variables ------------------------------ 706 INTEGER :: i 707 INTEGER :: j 708 INTEGER :: k 709 710 IF ( .NOT. this%para ) & 711 CALL ERROR("Ctqmc_allocateOpt : Ctqmc_setParameters never called ") 712 713 IF ( this%opt_analysis .EQ. 1 ) THEN 714 FREEIF(this%measCorrelation) 715 MALLOC(this%measCorrelation,(1:this%samples+1,1:3,1:this%flavors)) 716 this%measCorrelation = 0.d0 717 END IF 718 719 IF ( this%opt_order .GT. 0 ) THEN 720 FREEIF(this%measPerturbation) 721 MALLOC(this%measPerturbation,(1:this%opt_order,1:this%flavors)) 722 this%measPerturbation = 0.d0 723 END IF 724 725 IF ( this%opt_histo .GT. 0 ) THEN 726 FREEIF(this%occup_histo_time) 727 MALLOC(this%occup_histo_time,(1:this%flavors+1)) 728 this%occup_histo_time= 0.d0 729 FREEIF(this%occupconfig) 730 MALLOC(this%occupconfig,(1:2**this%flavors)) 731 this%occupconfig= 0.d0 732 FREEIF(this%suscep) 733 MALLOC(this%suscep,(1:3,1:this%samples)) 734 this%suscep= 0.d0 735 FREEIF(this%chi) 736 MALLOC(this%chi,(1:3,1:this%samples)) 737 this%chi= 0.d0 738 FREEIF(this%chicharge) 739 MALLOC(this%chicharge,(1:3,1:this%samples)) 740 this%chicharge= 0.d0 741 FREEIF(this%ntot) 742 MALLOC(this%ntot,(1:3)) 743 this%ntot= 0.d0 744 END IF 745 746 IF ( this%opt_noise .EQ. 1 ) THEN 747 IF ( ALLOCATED(this%measNoiseG) ) THEN 748 DO i=1,2 749 DO j = 1, this%flavors 750 DO k= 1, this%samples+1 751 CALL Vector_destroy(this%measNoiseG(k,j,i)) 752 END DO 753 END DO 754 END DO 755 DT_FREE(this%measNoiseG) 756 END IF 757 DT_MALLOC(this%measNoiseG,(1:this%samples+1,1:this%flavors,1:2)) 758 !DO i=1,2 759 DO j = 1, this%flavors 760 DO k= 1, this%samples+1 761 CALL Vector_init(this%measNoiseG(k,j,1),CTQMC_SLICE1) 762 END DO 763 END DO 764 DO j = 1, this%flavors 765 DO k= 1, this%samples+1 766 CALL Vector_init(this%measNoiseG(k,j,2),CTQMC_SLICE1*CTQMC_SLICE2+1) ! +1 pour etre remplacer ceil 767 END DO 768 END DO 769 !END DO 770 FREEIF(this%abNoiseG) 771 MALLOC(this%aBNoiseG,(1:2,1:this%samples+1,this%flavors)) 772 this%abNoiseG = 0.d0 773 END IF 774 775 IF (this%opt_spectra .GE. 1 ) THEN 776 FREEIF(this%density) 777 !MALLOC(this%density,(1:this%thermalization,1:this%flavors)) 778 i = CEILING(DBLE(this%thermalization+this%sweeps)/DBLE(this%measurements*this%opt_spectra)) 779 MALLOC(this%density,(1:this%flavors+1,1:i)) 780 this%density = 0.d0 781 END IF 782 !#endif 783 END SUBROUTINE Ctqmc_allocateOpt
ABINIT/m_Ctqmc/Ctqmc_clear [ Functions ]
NAME
Ctqmc_clear
FUNCTION
clear a ctqmc run
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
977 SUBROUTINE Ctqmc_clear(this) 978 979 !Arguments ------------------------------------ 980 TYPE(Ctqmc), INTENT(INOUT) :: this 981 !Local variables ------------------------------ 982 INTEGER :: i 983 INTEGER :: j 984 INTEGER :: k 985 986 this%measN(1,:) = 0.d0 987 this%measN(2,:) = 0.d0 988 !Do not set measN(3,:) to 0 to avoid erasing N between therm and ctqmc 989 this%measN(4,:) = 0.d0 990 this%measDE = 0.d0 991 ! this%seg_added = 0.d0 992 ! this%anti_added = 0.d0 993 ! this%seg_removed = 0.d0 994 ! this%anti_removed = 0.d0 995 ! this%seg_sign = 0.d0 996 ! this%anti_sign = 0.d0 997 this%stats(:) = 0.d0 998 this%swap = 0.d0 999 this%runTime = 0.d0 1000 this%modGlobalMove(2) = 0 1001 CALL Vector_clear(this%measNoise(1)) 1002 CALL Vector_clear(this%measNoise(2)) 1003 !#ifdef CTCtqmc_CHECK 1004 this%errorImpurity = 0.d0 1005 this%errorBath = 0.d0 1006 !#endif 1007 DO j = 1, this%flavors 1008 CALL GreenHyb_clear(this%Greens(j)) 1009 END DO 1010 !#ifdef CTCtqmc_ANALYSIS 1011 IF ( this%opt_analysis .EQ. 1 .AND. ALLOCATED(this%measCorrelation) ) & 1012 this%measCorrelation = 0.d0 1013 IF ( this%opt_order .GT. 0 .AND. ALLOCATED(this%measPerturbation) ) & 1014 this%measPerturbation = 0.d0 1015 IF ( this%opt_noise .EQ. 1 .AND. ALLOCATED(this%measNoiseG) ) THEN 1016 DO i=1,2 1017 DO j = 1, this%flavors 1018 DO k= 1, this%samples+1 1019 CALL Vector_clear(this%measNoiseG(k,j,i)) 1020 END DO 1021 END DO 1022 END DO 1023 END IF 1024 !#endif 1025 END SUBROUTINE Ctqmc_clear
ABINIT/m_Ctqmc/Ctqmc_computeF [ Functions ]
NAME
Ctqmc_computeF
FUNCTION
Compute the hybridization function
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc Gomega=G0 to compute F opt_fk=What is Gomega
OUTPUT
F=hybridization function
SIDE EFFECTS
NOTES
SOURCE
1158 SUBROUTINE Ctqmc_computeF(this, Gomega, F, opt_fk) 1159 1160 !Arguments ------------------------------------ 1161 TYPE(Ctqmc) , INTENT(INOUT) :: this 1162 COMPLEX(KIND=8), DIMENSION(:,:), INTENT(IN ) :: Gomega 1163 !INTEGER , INTENT(IN ) :: Wmax 1164 DOUBLE PRECISION, DIMENSION(:,:), INTENT(INOUT) :: F 1165 INTEGER , INTENT(IN ) :: opt_fk 1166 !Local variables ------------------------------ 1167 INTEGER :: flavors 1168 INTEGER :: samples 1169 INTEGER :: iflavor 1170 INTEGER :: iomega 1171 INTEGER :: itau 1172 DOUBLE PRECISION :: pi_invBeta 1173 DOUBLE PRECISION :: K 1174 DOUBLE PRECISION :: re 1175 DOUBLE PRECISION :: im 1176 COMPLEX(KIND=8), DIMENSION(:,:), ALLOCATABLE :: F_omega 1177 TYPE(GreenHyb) :: F_tmp 1178 1179 flavors = this%flavors 1180 1181 samples = this%samples 1182 pi_invBeta = ACOS(-1.d0) / this%beta 1183 this%Wmax=SIZE(Gomega,1) 1184 1185 IF ( this%have_MPI .EQV. .TRUE. ) THEN 1186 CALL GreenHyb_init(F_tmp,samples,this%beta,MY_COMM=this%MY_COMM) 1187 ELSE 1188 CALL GreenHyb_init(F_tmp,samples,this%beta) 1189 END IF 1190 ! K = this%mu 1191 1192 MALLOC(F_omega,(1:this%Wmax,1:flavors)) 1193 1194 !IF ( this%rank .EQ. 0 ) & 1195 !OPEN(UNIT=9876,FILE="K.dat",POSITION="APPEND") 1196 IF ( opt_fk .EQ. 0 ) THEN 1197 DO iflavor = 1, flavors 1198 DO iomega=1,this%Wmax 1199 re = REAL(Gomega(iomega,iflavor)) 1200 im = AIMAG(Gomega(iomega,iflavor)) 1201 F_omega(iomega,iflavor) = CMPLX(-re/(re*re+im*im),im/(re*re+im*im),8) 1202 END DO 1203 END DO 1204 !F_omega = CMPLX(-1.d0,0,8)/Gomega 1205 ELSE 1206 F_omega = Gomega 1207 END IF 1208 1209 DO iflavor = 1, flavors 1210 IF ( this%opt_levels .EQ. 1 ) THEN 1211 K = this%mu(iflavor) 1212 ELSE 1213 K = -REAL(F_omega(this%Wmax, iflavor)) 1214 ! this%mu = K 1215 this%mu(iflavor) = K 1216 END IF 1217 !IF ( this%rank .EQ. 0 ) & 1218 !WRITE(9876,'(I4,2E22.14)') iflavor, K, REAL(-F_omega(this%Wmax, iflavor)) 1219 !IF(this%rank .EQ.0) & 1220 !WRITE(this%ostream,*) "CTQMC K, this%mu = ",K,this%mu(iflavor) 1221 !WRITE(this%ostream,*) "CTQMC beta = ",this%beta 1222 IF ( opt_fk .EQ. 0 ) THEN 1223 DO iomega = 1, this%Wmax 1224 re = REAL(F_omega(iomega,iflavor)) 1225 im = AIMAG(F_omega(iomega,iflavor)) 1226 F_omega(iomega,iflavor) = CMPLX(re + K, im + (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, 8) 1227 !if(iflavor==1.and.this%rank==0) then 1228 !write(224,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(F_omega(iomega,iflavor)),imag(F_omega(iomega,iflavor)) 1229 !write(225,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(Gomega(iomega, iflavor)),imag(Gomega(iomega, iflavor)) 1230 !end if 1231 END DO 1232 ELSE 1233 DO iomega = 1, this%Wmax 1234 F_omega(iomega,iflavor) = F_omega(iomega,iflavor) & 1235 + CMPLX(K, 0.d0, 8) 1236 !if(iflavor==1.and.this%rank==0) then 1237 !write(224,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(F_omega(iomega,iflavor)),imag(F_omega(iomega,iflavor)) 1238 !write(225,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(Gomega(iomega, iflavor)),imag(Gomega(iomega, iflavor)) 1239 !end if 1240 END DO 1241 END IF 1242 K = REAL(CMPLX(0,(2.d0*DBLE(this%Wmax)-1.d0)*pi_invBeta,8)*F_omega(this%Wmax,iflavor)) 1243 CALL GreenHyb_setMuD1(this%Greens(iflavor),this%mu(iflavor),K) 1244 CALL GreenHyb_setOperW(F_tmp,F_omega(:,iflavor)) 1245 !CALL GreenHyb_backFourier(F_tmp,F_omega(:,iflavor)) 1246 CALL GreenHyb_backFourier(F_tmp) 1247 F(1:samples+1,iflavor) = (/ (-F_tmp%oper(samples+1-itau),itau=0,samples) /) 1248 END DO 1249 IF ( this%rank .EQ. 0 ) THEN 1250 open (unit=346,file='Hybridization.dat',status='unknown',form='formatted') 1251 DO iflavor = 1, flavors 1252 write(346,*) "#",iflavor 1253 do itau=1,this%samples+1 1254 write(346,*) itau,F(itau,iflavor) 1255 enddo 1256 write(346,*) 1257 END DO 1258 close(346) 1259 ENDIF 1260 FREE(F_omega) 1261 CALL GreenHyb_destroy(F_tmp) 1262 END SUBROUTINE Ctqmc_computeF
ABINIT/m_Ctqmc/Ctqmc_destroy [ Functions ]
NAME
Ctqmc_destroy
FUNCTION
destroy and deallocate all variables
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3577 SUBROUTINE Ctqmc_destroy(this) 3578 3579 !Arguments ------------------------------------ 3580 TYPE(Ctqmc), INTENT(INOUT) :: this 3581 !Local variables ------------------------------ 3582 INTEGER :: iflavor 3583 INTEGER :: flavors 3584 INTEGER :: i 3585 INTEGER :: j 3586 INTEGER :: k 3587 3588 if ( this%init .EQV. .FALSE. ) RETURN 3589 3590 flavors = this%flavors 3591 3592 CALL ImpurityOperator_destroy(this%Impurity) 3593 CALL BathOperator_destroy(this%Bath) 3594 CALL Vector_destroy(this%measNoise(1)) 3595 CALL Vector_destroy(this%measNoise(2)) 3596 3597 IF ( ALLOCATED(this%Greens) ) THEN 3598 DO iflavor = 1, flavors 3599 CALL GreenHyb_destroy(this%Greens(iflavor)) 3600 END DO 3601 DT_FREE( this%Greens ) 3602 END IF 3603 !#ifdef CTCtqmc_ANALYSIS 3604 FREEIF(this%measCorrelation) 3605 FREEIF(this%measPerturbation) 3606 IF ( this%opt_histo .GT. 0 ) THEN 3607 FREEIF(this%occup_histo_time) 3608 FREEIF(this%occupconfig) 3609 FREEIF(this%suscep) 3610 FREEIF(this%chi) 3611 FREEIF(this%chicharge) 3612 FREEIF(this%ntot) 3613 ENDIF 3614 FREEIF(this%measN) 3615 FREEIF(this%measDE) 3616 FREEIF(this%mu) 3617 FREEIF(this%abNoiseG) 3618 IF ( ALLOCATED(this%measNoiseG) ) THEN 3619 DO i=1,2 3620 DO j = 1, this%flavors 3621 DO k= 1, this%samples+1 3622 CALL Vector_destroy(this%measNoiseG(k,j,i)) 3623 END DO 3624 END DO 3625 END DO 3626 DT_FREE(this%measNoiseG) 3627 END IF 3628 FREEIF(this%density) 3629 !#endif 3630 this%ostream = 0 3631 this%istream = 0 3632 3633 this%sweeps = 0 3634 this%thermalization = 0 3635 this%flavors = 0 3636 this%samples = 0 3637 this%beta = 0.d0 3638 ! this%seg_added = 0.d0 3639 ! this%anti_added = 0.d0 3640 ! this%seg_removed = 0.d0 3641 ! this%anti_removed = 0.d0 3642 ! this%seg_sign = 0.d0 3643 ! this%anti_sign = 0.d0 3644 this%stats = 0.d0 3645 this%swap = 0.d0 3646 3647 3648 this%set = .FALSE. 3649 this%done = .FALSE. 3650 this%init = .FALSE. 3651 END SUBROUTINE Ctqmc_destroy
ABINIT/m_Ctqmc/Ctqmc_getD [ Functions ]
NAME
Ctqmc_getD
FUNCTION
get double occupation
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
D=full double occupation
SIDE EFFECTS
NOTES
SOURCE
2938 SUBROUTINE Ctqmc_getD(this, D) 2939 2940 !Arguments ------------------------------------ 2941 TYPE(Ctqmc) , INTENT(IN ) :: this 2942 DOUBLE PRECISION, DIMENSION(:,:), INTENT(OUT) :: D 2943 !Local variables ------------------------------ 2944 INTEGER :: iflavor1 2945 INTEGER :: iflavor2 2946 2947 IF ( SIZE(D,1) .LT. this%flavors .OR. SIZE(D,2) .LT. this%flavors ) & 2948 CALL ERROR("Ctqmc_getD : Dimensions of array D are too small") 2949 2950 D = 0.d0 2951 2952 DO iflavor1 = 1, this%flavors 2953 DO iflavor2 = 1, this%flavors 2954 D(iflavor2,iflavor1) = this%measDE(iflavor2,iflavor1) 2955 END DO 2956 D(iflavor1,iflavor1) = 0.d0 2957 END DO 2958 2959 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
E=interaction energy noise=noise on this value
SIDE EFFECTS
NOTES
SOURCE
2988 SUBROUTINE Ctqmc_getE(this,E,noise) 2989 2990 !Arguments ------------------------------------ 2991 TYPE(Ctqmc) , INTENT(IN ) :: this 2992 DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: E 2993 DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: Noise 2994 2995 IF ( PRESENT(E) ) & 2996 E = this%measDE(1,1) 2997 IF ( PRESENT(Noise) ) & 2998 Noise = SUM(this%Impurity%mat_U)/(this%flavors*(this%flavors-1)) & 2999 * this%a_Noise*(DBLE(this%sweeps)*DBLE(this%size))**this%b_Noise 3000 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
Gtau=green function in time Gw=green function in frequency
SIDE EFFECTS
NOTES
SOURCE
2857 SUBROUTINE Ctqmc_getGreen(this, Gtau, Gw) 2858 2859 !Arguments ------------------------------------ 2860 TYPE(Ctqmc) , INTENT(INOUT) :: this 2861 DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT(INOUT) :: Gtau 2862 COMPLEX(KIND=8), DIMENSION(:,:), OPTIONAL, INTENT(INOUT) :: Gw 2863 !Local variables ------------------------------ 2864 !INTEGER :: itime 2865 INTEGER :: iflavor1 2866 INTEGER :: iflavor2 2867 INTEGER :: iflavor3 2868 INTEGER :: flavors 2869 DOUBLE PRECISION :: u1 2870 DOUBLE PRECISION :: u2 2871 DOUBLE PRECISION :: Un 2872 DOUBLE PRECISION :: UUnn 2873 2874 flavors = this%flavors 2875 DO iflavor1 = 1, flavors 2876 u1 = 0.d0 2877 u2 = 0.d0 2878 DO iflavor2 = 1, flavors 2879 IF ( iflavor2 .EQ. iflavor1 ) CYCLE 2880 Un = this%Impurity%mat_U(iflavor2,iflavor1) * this%measN(1,iflavor2) 2881 u1 = u1 + Un 2882 u2 = u2 + Un*this%Impurity%mat_U(iflavor2,iflavor1) 2883 DO iflavor3 = 1, flavors 2884 IF ( iflavor3 .EQ. iflavor2 .OR. iflavor3 .EQ. iflavor1 ) CYCLE 2885 UUnn = (this%Impurity%mat_U(iflavor2,iflavor1)*this%Impurity%mat_U(iflavor3,iflavor1)) * this%measDE(iflavor2,iflavor3) 2886 u2 = u2 + UUnn 2887 END DO 2888 END DO 2889 2890 CALL GreenHyb_setMoments(this%Greens(iflavor1),u1,u2) 2891 IF ( PRESENT( Gtau ) ) THEN 2892 Gtau(1:this%samples,iflavor1) = this%Greens(iflavor1)%oper(1:this%samples) 2893 END IF 2894 !write(6,*) "present gw", present(gw) 2895 IF ( PRESENT( Gw ) ) THEN 2896 !write(6,*) "size gw",SIZE(Gw,DIM=2) ,flavors+1 2897 IF ( SIZE(Gw,DIM=2) .EQ. flavors+1 ) THEN 2898 CALL GreenHyb_forFourier(this%Greens(iflavor1), Gomega=Gw(:,iflavor1), omega=Gw(:,this%flavors+1)) 2899 !IF ( this%rank .EQ. 0 ) write(20,*) Gw(:,iflavor1) 2900 ELSE IF ( SIZE(Gw,DIM=2) .EQ. flavors ) THEN 2901 CALL GreenHyb_forFourier(this%Greens(iflavor1),Gomega=Gw(:,iflavor1)) 2902 ELSE 2903 CALL WARNALL("Ctqmc_getGreen : Gw is not valid ") 2904 CALL GreenHyb_forFourier(this%Greens(iflavor1),Wmax=this%Wmax) 2905 END IF 2906 ELSE 2907 CALL GreenHyb_forFourier(this%Greens(iflavor1),Wmax=this%Wmax) 2908 END IF 2909 END DO 2910 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-2024 ABINIT group (J. Bieder,F. Gendron) 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
2151 SUBROUTINE Ctqmc_getResult(this,Iatom,fname) 2152 2153 2154 #ifdef HAVE_MPI1 2155 include 'mpif.h' 2156 #endif 2157 !Arguments ------------------------------------ 2158 TYPE(Ctqmc) , INTENT(INOUT) :: this 2159 INTEGER, INTENT(IN ) :: Iatom 2160 character(len=fnlen), INTENT(IN) :: fname 2161 !Local variables ------------------------------ 2162 INTEGER :: iflavor 2163 INTEGER :: flavors 2164 INTEGER :: itau 2165 INTEGER :: endDensity 2166 DOUBLE PRECISION :: inv_flavors 2167 DOUBLE PRECISION :: a 2168 DOUBLE PRECISION :: b 2169 DOUBLE PRECISION :: r 2170 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: alpha 2171 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: beta 2172 DOUBLE PRECISION, DIMENSION(1:2) :: TabX 2173 DOUBLE PRECISION, DIMENSION(1:2) :: TabY 2174 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: freqs 2175 INTEGER, ALLOCATABLE, DIMENSION(:) :: counts 2176 INTEGER, ALLOCATABLE, DIMENSION(:) :: displs 2177 INTEGER, ALLOCATABLE, DIMENSION(:) :: occtot 2178 INTEGER, ALLOCATABLE, DIMENSION(:) :: spintot 2179 INTEGER, ALLOCATABLE, DIMENSION(:,:) :: occ 2180 INTEGER :: sp1,spinmax,spinmin,dspin,nelec,spin 2181 INTEGER :: spAll 2182 INTEGER :: last 2183 INTEGER :: n1 2184 INTEGER :: n2,n3,quotient,remainder,signe 2185 INTEGER :: debut 2186 ! INTEGER :: fin 2187 character(len=2) :: atomnb 2188 ! character(len=fnlen) :: tmpfile 2189 ! INTEGER :: unt 2190 #ifdef HAVE_MPI 2191 INTEGER :: ierr 2192 DOUBLE PRECISION, DIMENSION(1) :: rtime 2193 #endif 2194 DOUBLE PRECISION :: inv_size,sumh,sumtot 2195 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: buffer 2196 TYPE(FFTHyb) :: FFTmrka 2197 2198 #if defined HAVE_MPI && !defined HAVE_MPI2_INPLACE 2199 DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:) :: buffer1_out,freqs_buf 2200 DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:,:) :: buffer2_out 2201 #endif 2202 2203 IF ( .NOT. this%done ) & 2204 CALL ERROR("Ctqmc_getResult : Simulation not run ") 2205 2206 flavors = this%flavors 2207 inv_flavors = 1.d0 / DBLE(flavors) 2208 2209 2210 inv_size = 1.d0 / DBLE(this%size) 2211 sp1 = 0 2212 spAll = 0 2213 2214 !#ifdef CTCtqmc_CHECK 2215 IF ( this%opt_check .GT. 0 ) THEN 2216 this%errorImpurity = ImpurityOperator_getError(this%Impurity) * inv_flavors 2217 this%errorBath = BathOperator_getError (this%Bath ) * inv_flavors 2218 END IF 2219 !#endif 2220 2221 MALLOC(alpha,(1,1)) 2222 MALLOC(beta,(1,1)) 2223 MALLOC(buffer,(1,1)) 2224 IF ( this%opt_noise .EQ. 1) THEN 2225 FREEIF(alpha) 2226 MALLOC(alpha,(1:this%samples+1,1:flavors)) 2227 FREEIF(beta) 2228 MALLOC(beta,(1:this%samples+1,1:flavors)) 2229 END IF 2230 2231 IF ( this%have_MPI .EQV. .TRUE.) THEN 2232 sp1 = this%samples+1 2233 spALL = sp1 + flavors + 6 2234 2235 !#ifdef CTCtqmc_ANALYSIS 2236 IF ( this%opt_analysis .EQ. 1 ) & 2237 spAll = spAll + 3*sp1 2238 IF ( this%opt_order .GT. 0 ) & 2239 spAll = spAll + this%opt_order 2240 IF ( this%opt_noise .EQ. 1 ) & 2241 spAll = spAll + 2*(this%samples + 1) 2242 !#endif 2243 2244 FREEIF(buffer) 2245 MALLOC(buffer,(1:spAll,1:MAX(2,flavors))) 2246 END IF 2247 2248 ! this%seg_added = this%seg_added * inv_flavors 2249 ! this%seg_removed = this%seg_removed * inv_flavors 2250 ! this%seg_sign = this%seg_sign * inv_flavors 2251 ! this%anti_added = this%anti_added * inv_flavors 2252 ! this%anti_removed = this%anti_removed * inv_flavors 2253 ! this%anti_sign = this%anti_sign * inv_flavors 2254 this%stats(:) = this%stats(:) * inv_flavors 2255 2256 DO iflavor = 1, flavors 2257 CALL GreenHyb_measHybrid(this%Greens(iflavor), this%Bath%M(iflavor), this%Impurity%Particles(iflavor), .TRUE.) 2258 CALL GreenHyb_getHybrid(this%Greens(iflavor)) 2259 ! Accumule les dernieres mesure de N 2260 this%measN(1,iflavor) = this%measN(1,iflavor) + this%measN(3,iflavor)*this%measN(4,iflavor) 2261 this%measN(2,iflavor) = this%measN(2,iflavor) + this%measN(4,iflavor) 2262 ! Reduction 2263 this%measN(1,iflavor) = this%measN(1,iflavor) / ( this%measN(2,iflavor) * this%beta ) 2264 ! Correction 2265 CALL GreenHyb_setN(this%Greens(iflavor), this%measN(1,iflavor)) 2266 !#ifdef CTCtqmc_ANALYSIS 2267 IF ( this%opt_order .GT. 0 ) & 2268 this%measPerturbation(: ,iflavor) = this%measPerturbation(:,iflavor) & 2269 / SUM(this%measPerturbation(:,iflavor)) 2270 IF ( this%opt_analysis .EQ. 1 ) THEN 2271 this%measCorrelation (:,1,iflavor) = this%measCorrelation (:,1,iflavor) & 2272 / SUM(this%measCorrelation (:,1,iflavor)) & 2273 * this%inv_dt 2274 this%measCorrelation (:,2,iflavor) = this%measCorrelation (:,2,iflavor) & 2275 / SUM(this%measCorrelation (:,2,iflavor)) & 2276 * this%inv_dt 2277 this%measCorrelation (:,3,iflavor) = this%measCorrelation (:,3,iflavor) & 2278 / SUM(this%measCorrelation (:,3,iflavor)) & 2279 * this%inv_dt 2280 END IF 2281 !#endif 2282 IF ( this%opt_noise .EQ. 1 ) THEN 2283 TabX(1) = DBLE(this%modNoise2) 2284 TabX(2) = DBLE(this%modNoise1) 2285 DO itau = 1, this%samples+1 2286 this%measNoiseG(itau,iflavor,2)%vec = -this%measNoiseG(itau,iflavor,2)%vec*this%inv_dt & 2287 /(this%beta*DBLE(this%modNoise2)) 2288 this%measNoiseG(itau,iflavor,1)%vec = -this%measNoiseG(itau,iflavor,1)%vec*this%inv_dt & 2289 /(this%beta*DBLE(this%modNoise1)) 2290 n2 = this%measNoiseG(itau,iflavor,2)%tail 2291 TabY(1) = Stat_deviation(this%measNoiseG(itau,iflavor,2)%vec(1:n2))!*SQRT(n2/(n2-1)) 2292 n1 = this%measNoiseG(itau,iflavor,1)%tail 2293 TabY(2) = Stat_deviation(this%measNoiseG(itau,iflavor,1)%vec(1:n1))!*SQRT(n1/(n1-1)) 2294 CALL Stat_powerReg(TabX,SQRT(2.d0*LOG(2.d0))*TabY,alpha(itau,iflavor),beta(itau,iflavor),r) 2295 ! ecart type -> 60% 2296 ! largeur a mi-hauteur d'une gaussienne -> sqrt(2*ln(2))*sigma 2297 END DO 2298 END IF 2299 2300 IF ( this%have_MPI .EQV. .TRUE. ) THEN 2301 buffer(1:sp1, iflavor) = this%Greens(iflavor)%oper(1:sp1) 2302 END IF 2303 END DO 2304 last = sp1 2305 2306 2307 this%measDE(:,:) = this%measDE(:,:) * DBLE(this%measurements) /(DBLE(this%sweeps)*this%beta) 2308 IF ( this%opt_histo .GT. 0 ) THEN 2309 this%occup_histo_time(:) = this%occup_histo_time(:) / INT(this%sweeps/this%measurements) 2310 this%occupconfig(:) = this%occupconfig(:) / INT(this%sweeps/this%measurements) 2311 this%suscep(:,:) = this%suscep(:,:) / INT(this%sweeps/this%measurements) 2312 this%chi(:,:) = this%chi(:,:) / INT(this%sweeps/this%measurements) 2313 this%chicharge(:,:) = this%chicharge(:,:) / INT(this%sweeps/this%measurements) 2314 this%ntot(:) = this%ntot(:) / INT(this%sweeps/this%measurements) 2315 END IF 2316 ! write(6,*) "=== Histogram of occupations for complete simulation ====",INT(this%sweeps/this%measurements) 2317 ! sumh=0 2318 ! do n1=1,this%flavors+1 2319 ! write(6,'(i4,f10.4)') n1-1, this%occup_histo_time(n1) 2320 ! sumh=sumh+this%occup_histo_time(n1) 2321 ! enddo 2322 ! write(6,*) "=================================",sumh 2323 2324 n1 = this%measNoise(1)%tail 2325 n2 = this%measNoise(2)%tail 2326 2327 ! On utilise freqs comme tableau de regroupement 2328 ! Gather de Noise1 2329 IF ( this%have_MPI .EQV. .TRUE. ) THEN 2330 MALLOC(counts,(1:this%size)) 2331 MALLOC(displs,(1:this%size)) 2332 FREEIF(freqs) 2333 MALLOC(freqs,(1:this%size*n1)) 2334 freqs = 0.d0 2335 counts(:) = n1 2336 displs(:) = (/ ( iflavor*n1, iflavor=0, this%size-1 ) /) 2337 #ifdef HAVE_MPI 2338 #if defined HAVE_MPI2_INPLACE 2339 freqs(n1*this%rank+1:n1*(this%rank+1)) = this%measNoise(1)%vec(1:n1) 2340 CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_PRECISION, & 2341 freqs, counts, displs, & 2342 MPI_DOUBLE_PRECISION, this%MY_COMM, ierr) 2343 #else 2344 MALLOC(freqs_buf,(n1)) 2345 freqs_buf(1:n1)=this%measNoise(1)%vec(1:n1) 2346 CALL MPI_ALLGATHERV(freqs_buf, n1, MPI_DOUBLE_PRECISION, & 2347 freqs, counts, displs, & 2348 MPI_DOUBLE_PRECISION, this%MY_COMM, ierr) 2349 FREE(freqs_buf) 2350 #endif 2351 #endif 2352 n1 = this%size*n1 2353 CALL Vector_setSize(this%measNoise(1),n1) 2354 this%measNoise(1)%vec(1:n1) = freqs(:) 2355 ! Gather de Noise2 2356 FREE(freqs) 2357 MALLOC(freqs,(1:this%size*n2)) 2358 freqs = 0.d0 2359 counts(:) = n2 2360 displs(:) = (/ ( iflavor*n2, iflavor=0, this%size-1 ) /) 2361 #ifdef HAVE_MPI 2362 #if defined HAVE_MPI2_INPLACE 2363 freqs(n2*this%rank+1:n2*(this%rank+1)) = this%measNoise(2)%vec(1:n2) 2364 CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_PRECISION, & 2365 freqs, counts, displs, & 2366 MPI_DOUBLE_PRECISION, this%MY_COMM, ierr) 2367 #else 2368 MALLOC(freqs_buf,(n2)) 2369 freqs_buf(1:n2)=this%measNoise(2)%vec(1:n2) 2370 CALL MPI_ALLGATHERV(freqs_buf, n2, MPI_DOUBLE_PRECISION, & 2371 freqs, counts, displs, & 2372 MPI_DOUBLE_PRECISION, this%MY_COMM, ierr) 2373 FREE(freqs_buf) 2374 #endif 2375 #endif 2376 n2 = this%size*n2 2377 CALL Vector_setSize(this%measNoise(2),n2) 2378 this%measNoise(2)%vec(1:n2) = freqs(:) 2379 FREE(counts) 2380 FREE(displs) 2381 FREE(freqs) 2382 END IF 2383 !n1 = this%measNoise(1)%tail 2384 !n2 = this%measNoise(2)%tail 2385 2386 ! Transformation des paquets pour que ca fit a CTQMC_SLICE(1|2) 2387 IF ( n1 .GT. CTQMC_SLICE1 ) THEN 2388 itau = n1/CTQMC_SLICE1 2389 MALLOC(freqs,(1:n1/itau)) 2390 DO debut=1, n1/itau 2391 freqs(debut)=SUM(this%measNoise(1)%vec((debut-1)*itau+1:itau*debut)) 2392 END DO 2393 freqs(:) = freqs(:)/DBLE(itau) 2394 this%modNoise1 = this%modNoise1*itau 2395 n1 = n1/itau 2396 CALL Vector_setSize(this%measNoise(1),n1) 2397 this%measNoise(1)%vec(1:n1) = freqs(:) 2398 FREE(freqs) 2399 END IF 2400 IF ( n2 .GT. CTQMC_SLICE1*CTQMC_SLICE2 ) THEN 2401 itau = n2/(CTQMC_SLICE1*CTQMC_SLICE2) 2402 MALLOC(freqs,(1:n2/itau)) 2403 DO debut=1, n2/itau 2404 freqs(debut)=SUM(this%measNoise(2)%vec((debut-1)*itau+1:itau*debut)) 2405 END DO 2406 freqs(:) = freqs(:)/DBLE(itau) 2407 this%modNoise2 = this%modNoise2*itau 2408 n2 = n2/itau 2409 CALL Vector_setSize(this%measNoise(2),n2) 2410 this%measNoise(2)%vec(1:n2) = freqs(:) 2411 FREE(freqs) 2412 END IF 2413 ! On peut s'amuser avec nos valeur d'energies 2414 !MALLOC(TabX,(1:20)) 2415 !MALLOC(TabY,(1:20)) 2416 2417 TabX(1) = DBLE(this%modNoise2) 2418 TabX(2) = DBLE(this%modNoise1) 2419 2420 ! Il faut calculer pour chaque modulo 10 ecarts type sur les donnes acquises 2421 this%measNoise(1)%vec(1:n1) = this%measNoise(1)%vec(1:n1)/(this%beta*DBLE(this%modNoise1))*DBLE(this%measurements) 2422 this%measNoise(2)%vec(1:n2) = this%measNoise(2)%vec(1:n2)/(this%beta*DBLE(this%modNoise2))*DBLE(this%measurements) 2423 ! CALL Vector_print(this%measNoise(1),this%rank+70) 2424 ! CALL Vector_print(this%measNoise(2),this%rank+50) 2425 ! DO iflavor=1,10 2426 ! debut = (iflavor-1)*n2/10+1 2427 ! fin = iflavor*n2/10 2428 ! TabY(iflavor) = Stat_deviation(this%measNoise(2)%vec(debut:fin)) 2429 ! debut = (iflavor-1)*n1/10+1 2430 ! fin = iflavor*n1/10 2431 ! TabY(10+iflavor) = Stat_deviation(this%measNoise(1)%vec(debut:fin)) 2432 ! END DO 2433 !! TabY(1:n) = (this%measNoise(2)%vec(1:n) & 2434 !! ) 2435 !! !/(this%beta*DBLE(this%modNoise2))*DBLE(this%measurements) & 2436 !! !- this%measDE(1,1)) 2437 !! TabY(this%measNoise(2)%tail+1:n+this%measNoise(2)%tail) = (this%measNoise(1)%vec(1:n) & 2438 !! ) 2439 !! ! /(this%beta*DBLE(this%modNoise1))*DBLE(this%measurements) & 2440 !! ! - this%measDE(1,1)) 2441 ! IF ( this%rank .EQ. 0 ) THEN 2442 ! DO iflavor=1,20 2443 ! write(45,*) TabX(iflavor), TabY(iflavor) 2444 ! END DO 2445 ! END IF 2446 ! 2447 2448 2449 TabY(1) = Stat_deviation(this%measNoise(2)%vec(1:n2))!*SQRT(n2/(n2-1)) 2450 !! write(this%rank+10,*) TabX(2) 2451 !! write(this%rank+40,*) TabX(1) 2452 !! CALL Vector_print(this%measNoise(1),this%rank+10) 2453 !! CALL Vector_print(this%measNoise(2),this%rank+40) 2454 !! CLOSE(this%rank+10) 2455 !! CLOSE(this%rank+40) 2456 TabY(2) = Stat_deviation(this%measNoise(1)%vec(1:n1))!*SQRT(n1/(n1-1)) 2457 !! ! Ecart carre moyen ~ ecart type mais non biaise. Serait moins precis. Aucun 2458 ! impact sur la pente, juste sur l'ordonnee a l'origine. 2459 2460 CALL Stat_powerReg(TabX,SQRT(2.d0*LOG(2.d0))*TabY,a,b,r) 2461 ! FREE(TabX) 2462 ! FREE(TabY) 2463 ! ecart type -> 60% 2464 ! largeur a mi-hauteur d'une gaussienne -> sqrt(2*ln(2))*sigma 2465 2466 !this%measDE(1,1) = SUM(this%measNoise(1)%vec(1:this%measNoise(1)%tail))/(DBLE(this%measNoise(1)%tail*this%modNoise1)*this%beta) 2467 !this%measDE(2:flavors,1:flavors) = this%measDE(2:flavors,1:flavors) /(DBLE(this%sweeps)*this%beta) 2468 CALL ImpurityOperator_getErrorOverlap(this%Impurity,this%measDE) 2469 ! Add the difference between true calculation and quick calculation of the 2470 ! last sweep overlap to measDE(2,2) 2471 !this%measDE = this%measDE * DBLE(this%measurements) 2472 IF ( this%have_MPI .EQV. .TRUE. ) THEN 2473 IF ( this%opt_analysis .EQ. 1 ) THEN 2474 buffer(last+1:last+sp1,:) = this%measCorrelation(:,1,:) 2475 last = last + sp1 2476 buffer(last+1:last+sp1,:) = this%measCorrelation(:,2,:) 2477 last = last + sp1 2478 buffer(last+1:last+sp1,:) = this%measCorrelation(:,3,:) 2479 last = last + sp1 2480 END IF 2481 IF ( this%opt_order .GT. 0 ) THEN 2482 buffer(last+1:last+this%opt_order, :) = this%measPerturbation(:,:) 2483 last = last + this%opt_order 2484 END IF 2485 IF ( this%opt_noise .EQ. 1 ) THEN 2486 buffer(last+1:last+this%samples+1,:) = alpha(:,:) 2487 last = last + this%samples + 1 2488 buffer(last+1:last+this%samples+1,:) = beta(:,:) 2489 last = last + this%samples + 1 2490 END IF 2491 ! this%measDE(2,2) = a*EXP(b*LOG(DBLE(this%sweeps*this%size))) 2492 buffer(spall-(flavors+5):spAll-6,:) = this%measDE(:,:) 2493 ! buffer(spAll ,1) = this%seg_added 2494 ! buffer(spAll-1,1) = this%seg_removed 2495 ! buffer(spAll-2,1) = this%seg_sign 2496 ! buffer(spAll ,2) = this%anti_added 2497 ! buffer(spAll-1,2) = this%anti_removed 2498 ! buffer(spAll-2,2) = this%anti_sign 2499 buffer(spAll ,1) = this%stats(1) 2500 buffer(spAll-1,1) = this%stats(2) 2501 buffer(spAll-2,1) = this%stats(3) 2502 buffer(spAll ,2) = this%stats(4) 2503 buffer(spAll-1,2) = this%stats(5) 2504 buffer(spAll-2,2) = this%stats(6) 2505 buffer(spAll-3,1) = this%swap 2506 buffer(spAll-3,2) = DBLE(this%modGlobalMove(2)) 2507 buffer(spAll-4,1) = a 2508 buffer(spAll-4,2) = b 2509 !#ifdef CTCtqmc_CHECK 2510 buffer(spAll-5,1) = this%errorImpurity 2511 buffer(spAll-5,2) = this%errorBath 2512 !#endif 2513 2514 #ifdef HAVE_MPI 2515 CALL MPI_ALLREDUCE([this%runTime], rtime, 1, MPI_DOUBLE_PRECISION, MPI_MAX, this%MY_COMM, ierr) 2516 this%runTime=rtime(1) 2517 #if defined HAVE_MPI2_INPLACE 2518 CALL MPI_ALLREDUCE(MPI_IN_PLACE, buffer, spAll*flavors, & 2519 MPI_DOUBLE_PRECISION, MPI_SUM, this%MY_COMM, ierr) 2520 IF ( this%opt_histo .GT. 0 ) THEN 2521 CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%occup_histo_time, flavors+1, MPI_DOUBLE_PRECISION, MPI_SUM, & 2522 this%MY_COMM, ierr) 2523 CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%occupconfig, 2**flavors, MPI_DOUBLE_PRECISION, MPI_SUM, & 2524 this%MY_COMM, ierr) 2525 CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%suscep, 3*this%samples, MPI_DOUBLE_PRECISION, MPI_SUM, & 2526 this%MY_COMM, ierr) 2527 CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%chi, 3*this%samples, MPI_DOUBLE_PRECISION, MPI_SUM, & 2528 this%MY_COMM, ierr) 2529 CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%chicharge, 3*this%samples, MPI_DOUBLE_PRECISION, MPI_SUM, & 2530 this%MY_COMM, ierr) 2531 CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%ntot, 3, MPI_DOUBLE_PRECISION, MPI_SUM, & 2532 this%MY_COMM, ierr) 2533 END IF 2534 #else 2535 MALLOC(buffer2_out,(spAll,flavors)) 2536 CALL MPI_ALLREDUCE(buffer, buffer2_out, spAll*flavors, & 2537 MPI_DOUBLE_PRECISION, MPI_SUM, this%MY_COMM, ierr) 2538 buffer(1:spAll,1:flavors)=buffer2_out(1:spAll,1:flavors) 2539 FREE(buffer2_out) 2540 IF ( this%opt_histo .GT. 0 ) THEN 2541 MALLOC(buffer1_out,(flavors+1)) 2542 CALL MPI_ALLREDUCE(this%occup_histo_time, buffer1_out, flavors+1, MPI_DOUBLE_PRECISION, MPI_SUM, & 2543 this%MY_COMM, ierr) 2544 this%occup_histo_time(1:flavors+1)=buffer1_out(1:flavors+1) 2545 FREE(buffer1_out) 2546 END IF 2547 #endif 2548 #endif 2549 2550 2551 buffer = buffer * inv_size 2552 this%measDE(:,:) = buffer(spall-(flavors+5):spAll-6,:) 2553 ! this%seg_added = buffer(spAll ,1) 2554 ! this%seg_removed = buffer(spAll-1,1) 2555 ! this%seg_sign = buffer(spAll-2,1) 2556 ! this%anti_added = buffer(spAll ,2) 2557 ! this%anti_removed = buffer(spAll-1,2) 2558 ! this%anti_sign = buffer(spAll-2,2) 2559 this%stats(1) = buffer(spAll ,1) 2560 this%stats(2) = buffer(spAll-1,1) 2561 this%stats(3) = buffer(spAll-2,1) 2562 this%stats(4) = buffer(spAll ,2) 2563 this%stats(5) = buffer(spAll-1,2) 2564 this%stats(6) = buffer(spAll-2,2) 2565 this%swap = buffer(spAll-3,1) 2566 this%modGlobalMove(2) = NINT(buffer(spAll-3,2)) 2567 a = buffer(spAll-4,1) 2568 b = buffer(spAll-4,2) 2569 !!#ifdef CTCtqmc_CHECK 2570 this%errorImpurity= buffer(spAll-5,1) 2571 this%errorBath = buffer(spAll-5,2) 2572 !#endif 2573 2574 DO iflavor = 1, flavors 2575 this%Greens(iflavor)%oper = buffer(1:sp1 , iflavor) 2576 END DO 2577 last = sp1 2578 IF ( this%opt_analysis .EQ. 1 ) THEN 2579 this%measCorrelation(:,1,:) = buffer(last+1:last+sp1,:) 2580 last = last + sp1 2581 this%measCorrelation(:,2,:) = buffer(last+1:last+sp1,:) 2582 last = last + sp1 2583 this%measCorrelation(:,3,:) = buffer(last+1:last+sp1,:) 2584 last = last + sp1 2585 END IF 2586 IF ( this%opt_order .GT. 0 ) THEN 2587 this%measPerturbation(:,:) = buffer(last+1:last+this%opt_order, :) 2588 last = last + this%opt_order 2589 END IF 2590 IF ( this%opt_noise .EQ. 1 ) THEN 2591 alpha(:,:) = buffer(last+1:last+this%samples+1,:) 2592 last = last + this%samples + 1 2593 beta(:,:) = buffer(last+1:last+this%samples+1,:) 2594 last = last + this%samples + 1 2595 END IF 2596 END IF 2597 DO iflavor = 1, flavors 2598 ! complete DE this 2599 this%measDE(iflavor, iflavor+1:flavors) = this%measDE(iflavor+1:flavors,iflavor) 2600 END DO 2601 FREE(buffer) 2602 2603 IF ( this%opt_spectra .GE. 1 ) THEN 2604 endDensity = SIZE(this%density,2) 2605 IF ( this%density(1,endDensity) .EQ. -1.d0 ) & 2606 endDensity = endDensity - 1 2607 CALL FFTHyb_init(FFTmrka,endDensity,DBLE(this%thermalization)/DBLE(this%measurements*this%opt_spectra)) 2608 ! Not very Beauty 2609 MALLOC(freqs,(1:FFTmrka%size/2)) 2610 DO iflavor = 1, flavors 2611 ! mean value is removed to supress the continue composent 2612 CALL FFTHyb_setData(FFTmrka,this%density(iflavor,1:endDensity)/this%beta+this%Greens(iflavor)%oper(this%samples+1)) 2613 CALL FFTHyb_run(FFTmrka,1) 2614 CALL FFTHyb_getData(FFTmrka,endDensity,this%density(iflavor,:),freqs) 2615 END DO 2616 this%density(flavors+1,:) = -1.d0 2617 this%density(flavors+1,1:FFTmrka%size/2) = freqs 2618 CALL FFTHyb_destroy(FFTmrka) 2619 FREE(freqs) 2620 END IF 2621 2622 this%a_Noise = a 2623 this%b_Noise = b 2624 IF ( this%opt_noise .EQ. 1 ) THEN 2625 this%abNoiseG(1,:,:) = alpha 2626 this%abNoiseG(2,:,:) = beta 2627 END IF 2628 FREE(alpha) 2629 FREE(beta) 2630 2631 IF ( this%opt_histo .GT. 0 ) THEN 2632 write(this%ostream,*) "=== Histogram of occupations for complete simulation ====" 2633 ! write(6,*) "sumh over procs", sumh 2634 sumh=0 2635 do n1=1,this%flavors+1 2636 write(this%ostream,'(i4,f10.4)') n1-1, this%occup_histo_time(n1)/float(this%size) 2637 sumh=sumh+this%occup_histo_time(n1)/float(this%size) 2638 enddo 2639 write(this%ostream,'(a,f10.4)') " all" , sumh 2640 write(this%ostream,*) "=================================" 2641 2642 2643 MALLOC(occ,(2**this%flavors,1:flavors)) 2644 MALLOC(occtot,(2**this%flavors)) 2645 MALLOC(spintot,(2**this%flavors)) 2646 do n1=1,2**this%flavors 2647 ! Compute occupations of individual Orbitals 2648 n3=n1-1 2649 occtot(n1)=0 2650 spintot(n1)=0 2651 signe=1 2652 do n2=1,this%flavors 2653 remainder=modulo(n3,2) 2654 quotient=(n3-remainder)/2 2655 occ(n1,n2)=remainder 2656 n3=quotient 2657 occtot(n1)=occtot(n1)+occ(n1,n2) 2658 if(n2>=6) signe =-1 2659 !if(n2>=7) signe =0 2660 spintot(n1)=spintot(n1)+occ(n1,n2)*signe 2661 enddo 2662 this%occupconfig(n1)=this%occupconfig(n1)/float(this%size) 2663 enddo 2664 2665 write(this%ostream,*) "=== Histogram of occupations of configurations for complete simulation ====" 2666 sumh=0 2667 if(this%flavors==14) then 2668 do n1=1,2**this%flavors 2669 write(this%ostream,'(i4,14i2,f20.2)') n1, (occ(n1,n2),n2=1,14),this%occupconfig(n1) 2670 sumh=sumh+this%occupconfig(n1) 2671 enddo 2672 else if (this%flavors==10) then 2673 do n1=1,2**this%flavors 2674 write(this%ostream,'(i4,10i2,f20.2)') n1, (occ(n1,n2),n2=1,10),this%occupconfig(n1) 2675 sumh=sumh+this%occupconfig(n1) 2676 enddo 2677 end if 2678 write(this%ostream,'(a,f10.4)') " all" , sumh 2679 2680 sumtot=0 2681 do nelec=0,10 2682 spinmin=modulo(nelec,2) 2683 if(nelec<=5) spinmax=nelec 2684 if(nelec>=6) spinmax=10-nelec 2685 dspin=1 2686 write(this%ostream,*) "=== Histogram of occupations of configurations for total number of electrons",nelec 2687 do spin=spinmin,spinmax,dspin 2688 sumh=0 2689 do n1=1,2**this%flavors 2690 if(occtot(n1)==nelec.and.abs(spintot(n1))==spin) then 2691 sumh=sumh+this%occupconfig(n1) 2692 if(this%flavors==10) then 2693 write(this%ostream,'(i8,10i2,a,i2,i3,f10.4)') n1,(occ(n1,n2),n2=1,this%flavors)," ",occtot(n1),spintot(n1),& 2694 &this%occupconfig(n1) 2695 else if(this%flavors==14) then 2696 write(this%ostream,'(i8,14i2,a,i2,i3,f10.4)') n1,(occ(n1,n2),n2=1,this%flavors)," ",occtot(n1),spintot(n1),& 2697 &this%occupconfig(n1) 2698 end if 2699 endif 2700 enddo 2701 write(this%ostream,'(a,i4,a,i4,a,f10.4)') " === Sum of weights for",nelec," electrons and spin",spin," is ",sumh 2702 sumtot=sumtot+sumh 2703 enddo 2704 enddo 2705 write(this%ostream,'(a,f10.4)') "Full sum is",sumtot 2706 FREE(occ) 2707 FREE(occtot) 2708 FREE(spintot) 2709 2710 !============================== 2711 ! Print Susceptibilities 2712 !============================== 2713 write(atomnb, '("0",i1)') Iatom 2714 !Local Magnetic Susceptibility 2715 if(this%opt_histo .gt. 1) then 2716 !Scalar 2717 if(this%nspinor .eq. 1) then 2718 open(unit=735,file=trim(fname)//'_LocalSpinSuscept_atom_'//atomnb//'.dat',status='unknown',form='formatted') 2719 write(735,*) '#Tau Total t2g eg' 2720 do n1=1,this%samples 2721 this%suscep(:,n1)=this%suscep(:,n1)/float(this%size)/float(this%samples) 2722 write(735,'(1x,f14.8,2x,f12.8,2x,f12.8,2x,f12.8)') (n1-1)*this%beta/this%samples,(this%suscep(n2,n1),n2=1,3) 2723 enddo 2724 2725 else 2726 !SOC 2727 open (unit=735,file=trim(fname)//'_LocalMagnSuscept_atom_'//atomnb//'.dat',status='unknown',form='formatted') 2728 write(735,*) '#Tau Total Orbital Spin' 2729 do n1=1,this%samples 2730 this%chi(:,n1)=this%chi(:,n1)/float(this%size)/float(this%samples) 2731 write(735,'(1x,f14.8,2x,f12.8,2x,f12.8,2x,f12.8)') (n1-1)*this%beta/this%samples,(this%chi(n2,n1),n2=1,3) 2732 enddo 2733 endif 2734 close(unit=735) 2735 endif 2736 2737 !Local Charge Susceptibility 2738 if(this%opt_histo .gt. 2) then 2739 this%ntot(:)=this%ntot(:)/float(this%size)/float(this%samples) 2740 open (unit=735,file=trim(fname)//'_LocalChargeSuscept_atom_'//atomnb//'.dat',status='unknown',form='formatted') 2741 write(735,*) '#Tau Total <ntot> ' 2742 do n1=1,this%samples 2743 this%chicharge(1,n1)=(this%chicharge(1,n1)/float(this%size)/float(this%samples))-(this%ntot(1)*this%ntot(1)) 2744 !this%chicharge(2,n1)=(this%chicharge(2,n1)/float(this%size)/float(this%samples))-(this%ntot(2)*this%ntot(2)) 2745 !this%chicharge(3,n1)=(this%chicharge(3,n1)/float(this%size)/float(this%samples))-(this%ntot(3)*this%ntot(3)) 2746 !write(735 '(1x,f14.8,2x,f14.8,2x,f14.8,2x,f14.8,2x,f14.8)') (n1-1)*this%beta/this%samples,(this%chicharge(n2,n1),n2=1,3),this%ntot(1) 2747 write(735, '(1x,f14.8,2x,f14.8,2x,f14.8)') (n1-1)*this%beta/this%samples,(this%chicharge(1,n1)),this%ntot(1) 2748 enddo 2749 close(unit=735) 2750 endif 2751 2752 #if defined HAVE_FC_FLUSH 2753 call flush(735) 2754 #elif defined HAVE_FC_FLUSH_ 2755 call flush(735) 2756 #endif 2757 2758 2759 2760 2761 ENDIF 2762 2763 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=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
324 SUBROUTINE Ctqmc_init(this, ostream, istream, bFile, MY_COMM, iBuffer) 325 326 327 #ifdef HAVE_MPI1 328 include 'mpif.h' 329 #endif 330 !Arguments ------------------------------------ 331 TYPE(Ctqmc), INTENT(INOUT) :: this 332 INTEGER , INTENT(IN ) :: ostream 333 INTEGER , INTENT(IN ) :: istream 334 LOGICAL , INTENT(IN ) :: bFile 335 DOUBLE PRECISION, DIMENSION(1:10), OPTIONAL, INTENT(IN) :: iBuffer 336 INTEGER , OPTIONAL, INTENT(IN ) :: MY_COMM 337 !Local variables ------------------------------ 338 #ifdef HAVE_MPI 339 INTEGER :: ierr 340 #endif 341 INTEGER :: iflavor 342 #ifdef __GFORTRAN__ 343 ! INTEGER :: pid 344 ! CHARACTER(LEN=5) :: Cpid 345 ! 346 #endif 347 DOUBLE PRECISION, DIMENSION(1:10) :: buffer 348 349 this%ostream = ostream 350 this%istream = istream 351 352 ! --- RENICE --- 353 !#ifdef __GFORTRAN__ 354 ! pid = GetPid() 355 ! WRITE(Cpid,'(I5)') pid 356 ! CALL SYSTEM('renice +19 '//TRIM(ADJUSTL(Cpid))//' > /dev/null') 357 !#endif 358 !! --- RENICE --- 359 360 IF ( PRESENT(MY_COMM)) THEN 361 #ifdef HAVE_MPI 362 this%have_MPI = .TRUE. 363 this%MY_COMM = MY_COMM 364 CALL MPI_Comm_rank(this%MY_COMM, this%rank, ierr) 365 CALL MPI_Comm_size(this%MY_COMM, this%size, ierr) 366 #else 367 CALL WARN("Ctqmc_init : MPI is not used ") 368 this%have_MPI = .FALSE. 369 this%MY_COMM = -1 370 this%rank = 0 371 this%size = 1 372 #endif 373 ELSE 374 this%have_MPI = .FALSE. 375 this%MY_COMM = -1 376 this%rank = 0 377 this%size = 1 378 END IF 379 380 !IF ( this%rank .EQ. 0 ) THEN 381 ! WRITE(ostream,'(A20)') 'Job reniced with +19' 382 !CALL FLUSH(ostream) 383 !END IF 384 385 IF ( bFile .EQV. .TRUE. ) THEN 386 IF ( this%rank .EQ. 0 ) THEN 387 388 READ(istream,*) buffer(1) !iseed 389 READ(istream,*) buffer(2) !this%sweeps 390 READ(istream,*) buffer(3) !this%thermalization 391 READ(istream,*) buffer(4) !this%measurements 392 READ(istream,*) buffer(5) !this%flavors 393 READ(istream,*) buffer(6) !this%samples 394 READ(istream,*) buffer(7) !this%beta 395 READ(istream,*) buffer(8) !U 396 READ(istream,*) buffer(9) !iTech 397 READ(istream,*) buffer(10)!this%nspinor 398 !READ(istream,*) buffer(9) !Wmax 399 !#ifdef CTCtqmc_ANALYSIS 400 !READ(istream,*) buffer(10) !order 401 !#endif 402 END IF 403 404 #ifdef HAVE_MPI 405 IF ( this%have_MPI .EQV. .TRUE. ) & 406 CALL MPI_Bcast(buffer, 10, MPI_DOUBLE_PRECISION, 0, & 407 this%MY_COMM, ierr) 408 #endif 409 ELSE IF ( PRESENT(iBuffer) ) THEN 410 buffer(1:10) = iBuffer(1:10) 411 ELSE 412 CALL ERROR("Ctqmc_init : No input parameters ") 413 END IF 414 415 CALL Ctqmc_setParameters(this, buffer) 416 417 CALL Ctqmc_allocateAll(this) 418 419 DO iflavor = 1, this%flavors 420 CALL GreenHyb_init(this%Greens(iflavor),this%samples, this%beta, iTech=INT(buffer(9)),MY_COMM=this%MY_COMM) 421 END DO 422 423 424 ! this%seg_added = 0.d0 425 ! this%anti_added = 0.d0 426 ! this%seg_removed = 0.d0 427 ! this%anti_removed = 0.d0 428 ! this%seg_sign = 0.d0 429 ! this%anti_sign = 0.d0 430 this%stats(:) = 0.d0 431 this%swap = 0.d0 432 this%runTime = 0.d0 433 434 CALL Vector_init(this%measNoise(1),this%sweeps/this%modNoise1) 435 CALL Vector_init(this%measNoise(2),(this%sweeps/this%modNoise1+1)*CTQMC_SLICE2) 436 !CALL Vector_init(this%measNoise(3),101) 437 !CALL Vector_init(this%measNoise(4),101) 438 439 this%set = this%para .AND. this%inF 440 this%done = .FALSE. 441 this%init = .TRUE. 442 443 !#ifdef CTCtqmc_CHECK 444 this%errorImpurity = 0.d0 445 this%errorBath = 0.d0 446 !#endif 447 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc itotal=number of sweeps to perform : thermalization or sweeps ilatex=unit of file to write movie if so
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1504 SUBROUTINE Ctqmc_loop(this,itotal,ilatex) 1505 1506 !Arguments ------------------------------------ 1507 TYPE(Ctqmc), INTENT(INOUT) :: this 1508 INTEGER , INTENT(IN ) :: itotal 1509 INTEGER , INTENT(IN ) :: ilatex 1510 !Local variables ------------------------------ 1511 LOGICAL :: updated 1512 LOGICAL :: updated_seg 1513 LOGICAL, DIMENSION(:), ALLOCATABLE :: updated_swap 1514 1515 INTEGER :: flavors 1516 INTEGER :: measurements 1517 INTEGER :: modNoise1 1518 INTEGER :: modNoise2 1519 INTEGER :: modGlobalMove 1520 INTEGER :: sp1 1521 INTEGER :: itau 1522 INTEGER :: ind 1523 INTEGER :: endDensity 1524 INTEGER :: indDensity 1525 INTEGER :: swapUpdate1 1526 INTEGER :: swapUpdate2 1527 INTEGER :: old_percent 1528 INTEGER :: new_percent 1529 INTEGER :: ipercent 1530 INTEGER :: iflavor 1531 INTEGER :: isweep 1532 1533 DOUBLE PRECISION :: cpu_time1 1534 DOUBLE PRECISION :: cpu_time2 1535 DOUBLE PRECISION :: NRJ_old1 1536 DOUBLE PRECISION :: NRJ_old2 1537 DOUBLE PRECISION :: NRJ_new 1538 DOUBLE PRECISION :: total 1539 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_old1 1540 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_old2 1541 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_new 1542 1543 CALL CPU_TIME(cpu_time1) 1544 1545 flavors = this%flavors 1546 measurements = this%measurements 1547 modNoise1 = this%modNoise1 1548 modNoise2 = this%modNoise2 1549 modGlobalMove = this%modGlobalMove(1) 1550 sp1 = this%samples+1 1551 IF ( this%opt_histo .GT. 0 ) THEN 1552 this%occup_histo_time= 0.d0 1553 this%occupconfig= 0.d0 1554 this%suscep= 0.d0 1555 this%chi= 0.d0 1556 this%chicharge= 0.d0 1557 this%ntot = 0.d0 1558 END IF 1559 old_percent = 0 1560 1561 MALLOC(updated_swap,(1:flavors)) 1562 updated_swap(:) = .FALSE. 1563 1564 NRJ_old1 = 0.d0 1565 NRJ_old2 = 0.d0 1566 NRJ_new = 0.d0 1567 1568 MALLOC(gtmp_new,(1,1)) 1569 gtmp_new = 0.d0 1570 MALLOC(gtmp_old1,(1,1)) 1571 gtmp_old1 = 0.d0 1572 MALLOC(gtmp_old2,(1,1)) 1573 gtmp_old2 = 0.d0 1574 1575 endDensity = SIZE(this%density,2) 1576 1577 IF ( this%opt_noise .GT. 0 ) THEN 1578 FREEIF(gtmp_new) 1579 MALLOC(gtmp_new,(1:sp1,1:flavors)) 1580 FREEIF(gtmp_old1) 1581 MALLOC(gtmp_old1,(1:sp1,1:flavors)) 1582 FREEIF(gtmp_old2) 1583 MALLOC(gtmp_old2,(1:sp1,1:flavors)) 1584 END IF 1585 1586 IF ( this%rank .EQ. 0 ) THEN 1587 WRITE(this%ostream, '(1x,103A)') & 1588 "|----------------------------------------------------------------------------------------------------|" 1589 WRITE(this%ostream,'(1x,A)', ADVANCE="NO") "|" 1590 END IF 1591 1592 total = DBLE(itotal) 1593 1594 indDensity = 1 1595 DO isweep = 1, itotal 1596 DO iflavor = 1, flavors 1597 ImpurityOperator_QuickActivation(this%Impurity,iflavor) 1598 BathOperator_QuickActivation(this%Bath,iflavor) 1599 1600 CALL Ctqmc_tryAddRemove(this,updated_seg) 1601 updated = updated_seg .OR. updated_swap(iflavor) 1602 updated_swap(iflavor) = .FALSE. 1603 1604 CALL GreenHyb_measHybrid(this%Greens(iflavor), this%Bath%M(iflavor), this%Impurity%Particles(iflavor), updated) 1605 CALL Ctqmc_measN (this, iflavor, updated) 1606 IF ( this%opt_analysis .EQ. 1 ) & 1607 CALL Ctqmc_measCorrelation (this, iflavor) 1608 IF ( this%opt_order .GT. 0 ) & 1609 CALL Ctqmc_measPerturbation(this, iflavor) 1610 END DO 1611 1612 IF ( MOD(isweep,modGlobalMove) .EQ. 0 ) THEN 1613 CALL Ctqmc_trySwap(this,swapUpdate1, swapUpdate2) 1614 IF ( swapUpdate1 .NE. 0 .AND. swapUpdate2 .NE. 0 ) THEN 1615 updated_swap(swapUpdate1) = .TRUE. 1616 updated_swap(swapUpdate2) = .TRUE. 1617 END IF 1618 END IF 1619 1620 IF ( MOD(isweep,measurements) .EQ. 0 ) THEN 1621 CALL ImpurityOperator_measDE(this%Impurity,this%measDE) 1622 IF ( this%opt_spectra .GE. 1) THEN 1623 IF ( measurements*this%opt_spectra /= 0) THEN 1624 IF ( MOD(isweep,measurements*this%opt_spectra) .EQ. 0 ) THEN 1625 this%density(1:flavors,indDensity) = this%measN(3,1:flavors) 1626 indDensity = indDensity+1 1627 END IF 1628 END IF 1629 END IF 1630 END IF 1631 1632 IF ( MOD(isweep,measurements) .EQ. 0 ) THEN 1633 IF ( this%opt_histo .GT. 0 ) THEN 1634 CALL ImpurityOperator_occup_histo_time(this%Impurity,this%occup_histo_time,this%occupconfig,this%suscep,this%samples,& 1635 & this%chi,this%chicharge,this%ntot,this%opt_histo,this%nspinor) 1636 ENDIF 1637 ENDIF 1638 1639 IF ( MOD(isweep, modNoise1) .EQ. 0 ) THEN 1640 !modNext = isweep + modNoise2 1641 NRJ_new = (SUM(this%measDE(:,:))-this%measDE(1,1))*0.5d0 ! double occupation, avoid stat with 0 for U=J=0 1642 CALL Vector_pushBack(this%measNoise(1),NRJ_new - NRJ_old1) 1643 NRJ_old1 = NRJ_new 1644 1645 !! Try to limit accumulation error 1646 CALL ImpurityOperator_cleanOverlaps(this%Impurity) 1647 1648 IF ( this%opt_noise .EQ. 1 ) THEN 1649 DO iflavor = 1, flavors 1650 DO ind = 1, this%Greens(iflavor)%this%tail 1651 itau = this%Greens(iflavor)%this%listINT(ind) 1652 gtmp_new(itau,iflavor) = this%Greens(iflavor)%oper(itau) & 1653 +this%Greens(iflavor)%this%listDBLE(ind)*DBLE(this%Greens(iflavor)%factor) 1654 END DO 1655 DO itau = 1, sp1 1656 CALL Vector_pushBack(this%measNoiseG(itau,iflavor,1), gtmp_new(itau,iflavor) - gtmp_old1(itau,iflavor)) 1657 gtmp_old1(itau,iflavor) = gtmp_new(itau,iflavor) 1658 END DO 1659 END DO 1660 END IF 1661 END IF 1662 1663 IF ( MOD(isweep,modNoise2) .EQ. 0 ) THEN 1664 NRJ_new = (SUM(this%measDE(:,:))-this%measDE(1,1))*0.5d0 ! double occupation, avoid stat with 0 for U=J=0 1665 CALL Vector_pushBack(this%measNoise(2),NRJ_new - NRJ_old2) 1666 NRJ_old2 = NRJ_new 1667 IF ( this%opt_noise .EQ. 1 ) THEN 1668 DO iflavor = 1, flavors 1669 DO ind = 1, this%Greens(iflavor)%this%tail 1670 itau = this%Greens(iflavor)%this%listINT(ind) 1671 gtmp_new(itau,iflavor) = this%Greens(iflavor)%oper(itau) & 1672 +this%Greens(iflavor)%this%listDBLE(ind)*this%Greens(iflavor)%factor 1673 END DO 1674 DO itau = 1, sp1 1675 CALL Vector_pushBack(this%measNoiseG(itau,iflavor,2), gtmp_new(itau,iflavor) - gtmp_old2(itau,iflavor)) 1676 gtmp_old2(itau,iflavor) = gtmp_new(itau,iflavor) 1677 END DO 1678 END DO 1679 END IF 1680 1681 IF ( this%rank .EQ. 0 ) THEN 1682 new_percent = CEILING(DBLE(isweep)*100.d0/DBLE(itotal)) 1683 DO ipercent = old_percent+1, new_percent 1684 WRITE(this%ostream,'(A)',ADVANCE="NO") "-" 1685 END DO 1686 old_percent = new_percent 1687 END IF 1688 END IF 1689 1690 IF ( this%opt_movie .EQ. 1 ) THEN 1691 WRITE(ilatex,'(A11,I9)') "%iteration ", isweep 1692 CALL ImpurityOperator_printLatex(this%Impurity,ilatex,isweep) 1693 END IF 1694 1695 END DO 1696 1697 IF ( this%rank .EQ. 0 ) THEN 1698 DO ipercent = old_percent+1, 100 1699 WRITE(this%ostream,'(A)',ADVANCE="NO") "-" 1700 END DO 1701 WRITE(this%ostream,'(A)') "|" 1702 END IF 1703 1704 FREE(gtmp_new) 1705 FREE(gtmp_old1) 1706 FREE(gtmp_old2) 1707 FREE(updated_swap) 1708 1709 IF ( this%opt_spectra .GE. 1 .AND. itotal .EQ. this%sweeps ) THEN 1710 IF ( endDensity .NE. indDensity-1 ) THEN 1711 this%density(:,endDensity) = -1.d0 1712 END IF 1713 END IF 1714 1715 CALL CPU_TIME(cpu_time2) 1716 1717 this%runTime = (cpu_time2 - cpu_time1)*1.05d0 ! facteur arbitraire de correction 1718 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc iflavor=the flavor to measure
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2023 SUBROUTINE Ctqmc_measCorrelation(this, iflavor) 2024 2025 !Arguments ------------------------------------ 2026 TYPE(Ctqmc) , INTENT(INOUT) :: this 2027 !TYPE(ImpurityOperator), INTENT(IN ) :: impurity 2028 INTEGER , INTENT(IN ) :: iflavor 2029 !Local variables ------------------------------ 2030 INTEGER :: iCdag 2031 INTEGER :: iCdagBeta 2032 INTEGER :: iC 2033 INTEGER :: index 2034 INTEGER :: size 2035 DOUBLE PRECISION :: tC 2036 DOUBLE PRECISION :: tCdag 2037 !DOUBLE PRECISION :: time 2038 DOUBLE PRECISION :: inv_dt 2039 DOUBLE PRECISION :: beta 2040 2041 IF ( .NOT. this%set ) & 2042 CALL ERROR("Ctqmc_measCorrelation : QMC not set ") 2043 2044 size = this%impurity%particles(this%impurity%activeFlavor)%tail 2045 beta = this%beta 2046 2047 IF ( size .EQ. 0 ) RETURN 2048 2049 inv_dt = this%inv_dt 2050 2051 DO iCdag = 1, size ! first segments 2052 tCdag = this%impurity%particles(this%impurity%activeFlavor)%list(iCdag,Cdag_) 2053 tC = this%impurity%particles(this%impurity%activeFlavor)%list(iCdag,C_ ) 2054 index = INT( ( (tC - tCdag) * inv_dt ) + .5d0 ) + 1 2055 this%measCorrelation(index,1,iflavor) = this%measCorrelation(index,1,iflavor) + 1.d0 2056 MODCYCLE(iCdag+1,size,iCdagBeta) 2057 index = INT( ( ( & 2058 this%impurity%particles(this%impurity%activeFlavor)%list(iCdagBeta,Cdag_) - tC & 2059 + AINT(DBLE(iCdag)/DBLE(size))*beta & 2060 ) * inv_dt ) + .5d0 ) + 1 2061 IF ( index .LT. 1 .OR. index .GT. this%samples+1 ) THEN 2062 CALL WARN("Ctqmc_measCorrelation : bad index line 1095 ") 2063 ELSE 2064 this%measCorrelation(index,2,iflavor) = this%measCorrelation(index,2,iflavor) + 1.d0 2065 END IF 2066 ! DO iC = 1, size 2067 ! tC = impurity%particles(impurity%activeFlavor)%list(C_,iC) 2068 ! time = tC - tCdag 2069 ! IF ( time .LT. 0.d0 ) time = time + beta 2070 ! index = INT( ( time * inv_dt ) + .5d0 ) + 1 2071 ! this%measCorrelation(index,3,iflavor) = this%measCorrelation(index,3,iflavor) + 1.d0 2072 ! END DO 2073 DO iC = 1, size! this%Greens(iflavor)%index_old%tail 2074 this%measCorrelation(this%Greens(iflavor)%this%listINT(iC+(iCdag-1)*size),3,iflavor) = & 2075 this%measCorrelation(this%Greens(iflavor)%this%listINT(iC+(iCdag-1)*size),3,iflavor) + 1.d0 2076 END DO 2077 END DO 2078 2079 END SUBROUTINE Ctqmc_measCorrelation
ABINIT/m_Ctqmc/Ctqmc_measN [ Functions ]
NAME
Ctqmc_measN
FUNCTION
measure the number of electron
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc iflavor=which flavor to measure updated=something has changed since last time
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1974 SUBROUTINE Ctqmc_measN(this, iflavor, updated) 1975 1976 !Arguments ------------------------------------ 1977 TYPE(Ctqmc) , INTENT(INOUT) :: this 1978 !TYPE(ImpurityOperator), INTENT(IN ) :: impurity 1979 INTEGER , INTENT(IN ) :: iflavor 1980 LOGICAL , INTENT(IN ) :: updated 1981 1982 ! IF ( .NOT. this%set ) & 1983 ! CALL ERROR("Ctqmc_measN : QMC not set ") 1984 1985 1986 IF ( updated .EQV. .TRUE. ) THEN 1987 this%measN(1,iflavor) = this%measN(1,iflavor) + this%measN(3,iflavor)*this%measN(4,iflavor) 1988 this%measN(2,iflavor) = this%measN(2,iflavor) + this%measN(4,iflavor) 1989 this%measN(3,iflavor) = ImpurityOperator_measN(this%impurity) 1990 this%measN(4,iflavor) = 1.d0 1991 ELSE 1992 this%measN(4,iflavor) = this%measN(4,iflavor) + 1.d0 1993 END IF 1994 END SUBROUTINE Ctqmc_measN
ABINIT/m_Ctqmc/Ctqmc_measPerturbation [ Functions ]
NAME
Ctqmc_measPerturbation
FUNCTION
measure perturbation order
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc iflavor=the flavor to measure
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2107 SUBROUTINE Ctqmc_measPerturbation(this, iflavor) 2108 2109 !Arguments ------------------------------------ 2110 TYPE(Ctqmc) , INTENT(INOUT) :: this 2111 !TYPE(ImpurityOperator), INTENT(IN ) :: impurity 2112 INTEGER , INTENT(IN ) :: iflavor 2113 !Local variables ------------------------------ 2114 INTEGER :: index 2115 2116 IF ( .NOT. this%set ) & 2117 CALL ERROR("Ctqmc_measiPerturbation : QMC not set ") 2118 2119 index = this%impurity%particles(this%impurity%activeFlavor)%tail + 1 2120 IF ( index .LE. this%opt_order ) & 2121 this%measPerturbation(index,iflavor) = this%measPerturbation(index,iflavor) + 1.d0 2122 2123 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3027 SUBROUTINE Ctqmc_printAll(this) 3028 3029 !Arguments ------------------------------------ 3030 TYPE(Ctqmc), INTENT(INOUT) :: this 3031 3032 IF ( .NOT. this%done ) & 3033 CALL WARNALL("Ctqmc_printAll : Simulation not run ") 3034 3035 CALL Ctqmc_printQMC(this) 3036 3037 CALL Ctqmc_printGreen(this) 3038 3039 CALL Ctqmc_printD(this) 3040 3041 ! CALL Ctqmc_printE(this) 3042 3043 !#ifdef CTCtqmc_ANALYSIS 3044 CALL Ctqmc_printPerturbation(this) 3045 3046 CALL Ctqmc_printCorrelation(this) 3047 !#endif 3048 3049 CALL Ctqmc_printSpectra(this) 3050 3051 END SUBROUTINE Ctqmc_printAll
ABINIT/m_Ctqmc/Ctqmc_printCorrelation [ Functions ]
NAME
Ctqmc_printCorrelation
FUNCTION
print correlation fonctions
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3435 SUBROUTINE Ctqmc_printCorrelation(this, oFileIn) 3436 3437 !Arguments ------------------------------------ 3438 TYPE(Ctqmc) , INTENT(IN) :: this 3439 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3440 !Local variables ------------------------------ 3441 INTEGER :: oFile 3442 INTEGER :: itime 3443 INTEGER :: sp1 3444 INTEGER :: iflavor 3445 INTEGER :: i 3446 INTEGER :: flavors 3447 CHARACTER(LEN=2) :: a 3448 CHARACTER(LEN=50) :: string 3449 DOUBLE PRECISION :: dt 3450 3451 !IF ( this%rank .NE. MOD(5,this%size)) RETURN 3452 IF ( this%rank .NE. MOD(this%size+5,this%size)) RETURN 3453 IF ( this%opt_analysis .NE. 1 ) RETURN 3454 3455 oFile = 44 3456 IF ( PRESENT(oFileIn) ) THEN 3457 oFile = oFileIn 3458 ELSE 3459 OPEN(UNIT=oFile, FILE="Correlation.dat") 3460 END IF 3461 3462 sp1 = this%samples 3463 dt = this%beta / sp1 3464 sp1 = sp1 + 1 3465 flavors = this%flavors 3466 3467 i = 3*flavors + 1 3468 WRITE(a,'(I2)') i 3469 WRITE(oFile,*) "# time (/ (segement, antiseg, correl), i=1, flavor/)" 3470 string = '(1x,'//TRIM(ADJUSTL(a))//'F19.15)' 3471 DO itime = 1, sp1 3472 WRITE(oFile,string) DBLE(itime-1)*dt, & 3473 (/ ( & 3474 (/ ( this%measCorrelation(itime, i, iflavor), i=1,3) /) & 3475 , iflavor=1, flavors) /) 3476 END DO 3477 3478 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3479 3480 END SUBROUTINE Ctqmc_printCorrelation
ABINIT/m_Ctqmc/Ctqmc_printD [ Functions ]
NAME
Ctqmc_printD
FUNCTION
print individual double occupancy
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3259 SUBROUTINE Ctqmc_printD(this,oFileIn) 3260 3261 !Arguments ------------------------------------ 3262 TYPE(Ctqmc) , INTENT(IN) :: this 3263 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3264 !Local variables ------------------------------ 3265 INTEGER :: oFile 3266 INTEGER :: iflavor1 3267 INTEGER :: iflavor2 3268 3269 !IF ( this%rank .NE. MOD(2,this%size)) RETURN 3270 IF ( this%rank .NE. MOD(this%size+2,this%size)) RETURN 3271 3272 oFile = 41 3273 IF ( PRESENT(oFileIn) ) THEN 3274 oFile = oFileIn 3275 ELSE 3276 OPEN(UNIT=oFile, FILE="D.dat") 3277 END IF 3278 3279 DO iflavor1 = 1, this%flavors 3280 DO iflavor2 = iflavor1+1, this%flavors 3281 WRITE(oFile,'(1x,A8,I4,A1,I4,A3,ES21.14)') "Orbitals", iflavor1, "-", iflavor2, " : ", this%measDE(iflavor2,iflavor1) 3282 END DO 3283 END DO 3284 3285 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3286 3287 END SUBROUTINE Ctqmc_printD
ABINIT/m_Ctqmc/Ctqmc_printE [ Functions ]
NAME
Ctqmc_printE
FUNCTION
print energy and noise
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3315 SUBROUTINE Ctqmc_printE(this,oFileIn) 3316 3317 !Arguments ------------------------------------ 3318 TYPE(Ctqmc) , INTENT(IN) :: this 3319 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3320 !Local variables ------------------------------ 3321 INTEGER :: oFile 3322 DOUBLE PRECISION :: E 3323 DOUBLE PRECISION :: Noise 3324 3325 !IF ( this%rank .NE. MOD(3,this%size)) RETURN 3326 IF ( this%rank .NE. MOD(this%size+3,this%size)) RETURN 3327 3328 oFile = 42 3329 IF ( PRESENT(oFileIn) ) THEN 3330 oFile = oFileIn 3331 ELSE 3332 OPEN(UNIT=oFile, FILE="BetaENoise.dat") 3333 END IF 3334 3335 CALL Ctqmc_getE(this,E,Noise) 3336 3337 WRITE(oFile,'(1x,F5.2,A2,ES21.14,A2,ES21.14)') this%beta, " ", E, " ", Noise 3338 3339 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3340 3341 END SUBROUTINE Ctqmc_printE
ABINIT/m_Ctqmc/Ctqmc_printGreen [ Functions ]
NAME
Ctqmc_printGreen
FUNCTION
print green functions
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3180 SUBROUTINE Ctqmc_printGreen(this, oFileIn) 3181 3182 !Arguments ------------------------------------ 3183 TYPE(Ctqmc) , INTENT(IN) :: this 3184 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3185 !Local variables ------------------------------ 3186 INTEGER :: oFile 3187 INTEGER :: itime 3188 INTEGER :: sp1 3189 INTEGER :: iflavor 3190 INTEGER :: flavors 3191 CHARACTER(LEN=4) :: cflavors 3192 CHARACTER(LEN=50) :: string 3193 DOUBLE PRECISION :: dt 3194 DOUBLE PRECISION :: sweeps 3195 3196 !IF ( this%rank .NE. MOD(1,this%size)) RETURN 3197 IF ( this%rank .NE. MOD(this%size+1,this%size)) RETURN 3198 3199 oFile = 40 3200 IF ( PRESENT(oFileIn) ) THEN 3201 oFile = oFileIn 3202 ELSE 3203 OPEN(UNIT=oFile, FILE="Gtau.dat") 3204 END IF 3205 3206 sp1 = this%samples 3207 dt = this%beta / DBLE(sp1) 3208 sp1 = sp1 + 1 3209 flavors = this%flavors 3210 sweeps = DBLE(this%sweeps)*DBLE(this%size) 3211 3212 IF ( this%opt_noise .EQ. 1) THEN 3213 WRITE(cflavors,'(I4)') 2*flavors+1 3214 string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)' 3215 DO itime = 1, sp1 3216 WRITE(oFile,string) DBLE(itime-1)*dt, & 3217 (/ (this%Greens(iflavor)%oper(itime), iflavor=1, flavors) /), & 3218 (/ (this%abNoiseG(1,itime,iflavor)*(sweeps)**this%abNoiseG(2,itime,iflavor), iflavor=1, flavors) /) 3219 END DO 3220 ELSE 3221 WRITE(cflavors,'(I4)') flavors+1 3222 string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)' 3223 DO itime = 1, sp1 3224 WRITE(oFile,string) DBLE(itime-1)*dt, & 3225 (/ (this%Greens(iflavor)%oper(itime), iflavor=1, flavors) /) 3226 END DO 3227 END IF 3228 3229 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3230 3231 END SUBROUTINE Ctqmc_printGreen
ABINIT/m_Ctqmc/Ctqmc_printPerturbation [ Functions ]
NAME
Ctqmc_printPerturbation
FUNCTION
print perturbation order
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc oFileIn=file stream
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
3371 SUBROUTINE Ctqmc_printPerturbation(this, oFileIn) 3372 3373 !Arguments ------------------------------------ 3374 TYPE(Ctqmc) , INTENT(IN) :: this 3375 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3376 !Local variables------------------------------- 3377 INTEGER :: oFile 3378 INTEGER :: iorder 3379 INTEGER :: order 3380 INTEGER :: iflavor 3381 INTEGER :: flavors 3382 CHARACTER(LEN=2) :: a 3383 CHARACTER(LEN=50) :: string 3384 3385 !IF ( this%rank .NE. MOD(4,this%size)) RETURN 3386 IF ( this%rank .NE. MOD(this%size+4,this%size)) RETURN 3387 IF ( this%opt_order .LE. 0 ) RETURN 3388 3389 oFile = 43 3390 IF ( PRESENT(oFileIn) ) THEN 3391 oFile = oFileIn 3392 ELSE 3393 OPEN(UNIT=oFile, FILE="Perturbation.dat") 3394 END IF 3395 3396 order = this%opt_order 3397 flavors = this%flavors 3398 3399 WRITE(a,'(I2)') flavors 3400 string = '(I5,'//TRIM(ADJUSTL(a))//'F19.15)' 3401 DO iorder = 1, order 3402 WRITE(oFile,string) iorder-1, & 3403 (/ (this%measPerturbation(iorder, iflavor), iflavor=1, flavors) /) 3404 END DO 3405 3406 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3407 END SUBROUTINE Ctqmc_printPerturbation
ABINIT/m_Ctqmc/Ctqmc_printQMC [ Functions ]
NAME
Ctqmc_printQMC
FUNCTION
print ctqmc statistics
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3078 SUBROUTINE Ctqmc_printQMC(this) 3079 3080 !Arguments ------------------------------------ 3081 TYPE(Ctqmc), INTENT(INOUT) :: this 3082 !Local variables ------------------------------ 3083 INTEGER :: ostream 3084 INTEGER :: iflavor 3085 DOUBLE PRECISION :: sweeps 3086 DOUBLE PRECISION :: invSweeps 3087 CHARACTER(LEN=2) :: a 3088 CHARACTER(LEN=15) :: string 3089 3090 !IF ( this%rank .NE. 0) RETURN 3091 IF ( this%rank .NE. MOD(this%size,this%size)) RETURN 3092 3093 ostream = this%ostream 3094 sweeps = DBLE(this%sweeps) 3095 invSweeps = 1.d0/sweeps 3096 3097 WRITE(ostream,'(1x,F13.0,A11,F10.2,A12,I5,A5)') sweeps*DBLE(this%size), " sweeps in ", this%runTime, & 3098 " seconds on ", this%size, " CPUs" 3099 WRITE(ostream,'(A28,F6.2)') "Segments added [%] : ", this%stats(CTQMC_SEGME+CTQMC_ADDED)*invSweeps*100.d0 3100 WRITE(ostream,'(A28,F6.2)') "Segments removed [%] : ", this%stats(CTQMC_SEGME+CTQMC_REMOV)*invSweeps*100.d0 3101 WRITE(ostream,'(A28,F6.2)') "Segments sign [%] : ", this%stats(CTQMC_SEGME+CTQMC_DETSI)*invSweeps*100.d0 3102 WRITE(ostream,'(A28,F6.2)') "Anti-segments added [%] : ", this%stats(CTQMC_ANTIS+CTQMC_ADDED)*invSweeps*100.d0 3103 WRITE(ostream,'(A28,F6.2)') "Anti-segments removed [%] : ", this%stats(CTQMC_ANTIS+CTQMC_REMOV)*invSweeps*100.d0 3104 WRITE(ostream,'(A28,F6.2)') "Anti-segments sign [%] : ", this%stats(CTQMC_ANTIS+CTQMC_DETSI)*invSweeps*100.d0 3105 IF ( this%modGlobalMove(1) .LT. this%sweeps + 1 ) THEN 3106 WRITE(ostream,'(A28,F6.2)') "Global Move [%] : ", this%swap *invSweeps*100.d0*this%modGlobalMove(1) 3107 WRITE(ostream,'(A28,F6.2)') "Global Move Reduced [%] : ", this%swap / DBLE(this%modGlobalMove(2))*100.d0 3108 END IF 3109 !#ifdef CTCtqmc_CHECK 3110 IF ( this%opt_check .EQ. 1 .OR. this%opt_check .EQ. 3 ) & 3111 WRITE(ostream,'(A28,E22.14)') "Impurity test [%] : ", this%errorImpurity*100.d0 3112 IF ( this%opt_check .GE. 2 ) & 3113 WRITE(ostream,'(A28,E22.14)') "Bath test [%] : ", this%errorBath *100.d0 3114 !#endif 3115 WRITE(ostream,'(A28,ES22.14,A5,ES21.14)') "<Epot> [U] : ", this%measDE(1,1), " +/- ",& 3116 !#ifdef HAVE_MPI 3117 SUM(this%Impurity%mat_U)/(this%flavors*(this%flavors-1)) * this%a_Noise*(sweeps*DBLE(this%size))**this%b_Noise 3118 !#else 3119 ! this%a_Noise*(sweeps)**this%b_Noise 3120 !#endif 3121 WRITE(ostream,'(A28,F8.4,A3,F7.4)') "Noise [/U] : ", this%a_Noise, " x^", this%b_Noise 3122 WRITE(ostream,'(A28,E10.2)') "Niquist puls. [/beta] : ", ACOS(-1.d0)*this%inv_dt 3123 WRITE(ostream,'(A28,E22.14)') "Max Acc. Epot Error [U] : ", this%measDE(2,2)/(this%beta*this%modNoise1*2.d0)*sweeps 3124 3125 !WRITE(ostream,'(A28,F7.4,A3,F7.4,A4,E20.14)') "Noise [G(tau)] : ", this%a_Noise(2), "x^", this%b_Noise(2), " -> ", & 3126 !this%a_Noise(2)*(sweeps*DBLE(this%size))**this%b_Noise(2) 3127 IF ( this%opt_order .GT. 0 ) THEN 3128 WRITE(a,'(I2)') this%flavors 3129 string = '(A28,'//TRIM(ADJUSTL(a))//'(1x,I3))' 3130 WRITE(ostream,string) "Perturbation orders : ", & 3131 (/ (MAXLOC(this%measPerturbation(:, iflavor))-1, iflavor=1, this%flavors) /) 3132 END IF 3133 !CALL FLUSH(this%ostream) 3134 IF ( ABS(((this%stats(CTQMC_SEGME+CTQMC_ADDED) *invSweeps*100.d0) / & 3135 (this%stats(CTQMC_SEGME+CTQMC_REMOV) *invSweeps*100.d0) - 1.d0)) .GE. 0.02d0 & 3136 .OR. ABS(((this%stats(CTQMC_ANTIS+CTQMC_ADDED)*invSweeps*100.d0) / & 3137 (this%stats(CTQMC_ANTIS+CTQMC_REMOV)*invSweeps*100.d0) - 1.d0)) .GE. 0.02d0 ) & 3138 THEN 3139 CALL WARNALL("Ctqmc_printQMC : bad statistic according to moves. Increase sweeps") 3140 END IF 3141 ! Check sign problem for diagonal hybridization. 3142 IF ( (this%stats(CTQMC_SEGME+CTQMC_DETSI) + this%stats(CTQMC_ANTIS+CTQMC_DETSI)) .GT. 1.d-10 ) THEN 3143 CALL WARNALL("Ctqmc_printQMC : at least one negative sign occured. There might be a bug in the CT-QMC") 3144 END IF 3145 3146 IF ( ABS(this%b_Noise+0.5)/0.5d0 .GE. 0.05d0 ) & 3147 CALL WARNALL("Ctqmc_printQMC : bad statistic according to Noise. Increase sweeps") 3148 ! IF ( ISNAN(this%a_Noise) .OR. ISNAN(this%a_Noise) ) & 3149 ! CALL WARNALL("Ctqmc_printQMC : NaN appeared. Increase sweeps ") 3150 3151 3152 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc oFileIn=file stream
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
3509 SUBROUTINE Ctqmc_printSpectra(this, oFileIn) 3510 3511 !Arguments ------------------------------------ 3512 TYPE(Ctqmc) , INTENT(IN) :: this 3513 INTEGER , OPTIONAL, INTENT(IN) :: oFileIn 3514 !Local variables ------------------------------ 3515 INTEGER :: oFile 3516 INTEGER :: flavors 3517 INTEGER :: indDensity 3518 INTEGER :: endDensity 3519 CHARACTER(LEN=4) :: a 3520 CHARACTER(LEN=16) :: formatSpectra 3521 3522 !IF ( this%rank .NE. MOD(6,this%size)) RETURN 3523 IF ( this%opt_spectra .LT. 1 ) RETURN 3524 3525 oFile = 45+this%rank 3526 a ="0000" 3527 WRITE(a,'(I4)') this%rank 3528 IF ( PRESENT(oFileIn) ) THEN 3529 oFile = oFileIn 3530 ELSE 3531 OPEN(UNIT=oFile, FILE="Markov_"//TRIM(ADJUSTL(a))//".dat") 3532 END IF 3533 3534 flavors = this%flavors 3535 WRITE(a,'(I4)') flavors+1 3536 formatSpectra ='(1x,'//TRIM(ADJUSTL(a))//'ES22.14)' 3537 WRITE(oFile,*) "# freq[/hermalization] FFT" 3538 3539 endDensity = SIZE(this%density,2) 3540 DO WHILE ( this%density(flavors+1,endDensity) .EQ. -1 ) 3541 endDensity = endDensity -1 3542 END DO 3543 3544 DO indDensity = 1, endDensity 3545 WRITE(oFile,formatSpectra) this%density(flavors+1,indDensity), this%density(1:flavors,indDensity) 3546 END DO 3547 3548 IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile) 3549 3550 END SUBROUTINE Ctqmc_printSpectra
ABINIT/m_Ctqmc/Ctqmc_reset [ Functions ]
NAME
Ctqmc_reset
FUNCTION
reset a ctqmc simulation
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1052 SUBROUTINE Ctqmc_reset(this) 1053 1054 !Arguments ------------------------------------ 1055 TYPE(Ctqmc), INTENT(INOUT) :: this 1056 !Local variables ------------------------------ 1057 INTEGER :: iflavor 1058 DOUBLE PRECISION :: sweeps 1059 1060 DO iflavor = 1, this%flavors 1061 CALL GreenHyb_reset(this%Greens(iflavor)) 1062 END DO 1063 CALL Ctqmc_clear(this) 1064 CALL ImpurityOperator_reset(this%Impurity) 1065 CALL BathOperator_reset (this%Bath) 1066 this%measN(3,:) = 0.d0 1067 !complete restart -> measN=0 1068 this%done = .FALSE. 1069 this%set = .FALSE. 1070 this%inF = .FALSE. 1071 this%opt_movie = 0 1072 this%opt_analysis = 0 1073 this%opt_order = 0 1074 this%opt_check = 0 1075 this%opt_noise = 0 1076 this%opt_spectra = 0 1077 this%opt_levels = 0 1078 sweeps = DBLE(this%sweeps)*DBLE(this%size) 1079 CALL Ctqmc_setSweeps(this, sweeps) 1080 !#ifdef HAVE_MPI 1081 ! CALL MPI_BARRIER(this%MY_COMM,iflavor) 1082 ! IF ( this%rank .EQ. 0 ) & 1083 !#endif 1084 ! WRITE(this%ostream,'(A9)') "QMC reset" 1085 ! CALL FLUSH(this%ostream) 1086 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=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
1345 SUBROUTINE Ctqmc_run(this,opt_order,opt_histo,opt_movie,opt_analysis,opt_check,opt_noise,opt_spectra,opt_gMove) 1346 1347 1348 #ifdef HAVE_MPI1 1349 include 'mpif.h' 1350 #endif 1351 !Arguments ------------------------------------ 1352 TYPE(Ctqmc), INTENT(INOUT) :: this 1353 INTEGER, OPTIONAL, INTENT(IN ) :: opt_order 1354 INTEGER, OPTIONAL, INTENT(IN ) :: opt_histo 1355 INTEGER, OPTIONAL, INTENT(IN ) :: opt_movie 1356 INTEGER, OPTIONAL, INTENT(IN ) :: opt_analysis 1357 INTEGER, OPTIONAL, INTENT(IN ) :: opt_check 1358 INTEGER, OPTIONAL, INTENT(IN ) :: opt_noise 1359 INTEGER, OPTIONAL, INTENT(IN ) :: opt_spectra 1360 INTEGER, OPTIONAL, INTENT(IN ) :: opt_gMove 1361 !Local variables ------------------------------ 1362 #ifdef HAVE_MPI 1363 INTEGER :: ierr 1364 #endif 1365 !#ifdef CTCtqmc_MOVIE 1366 INTEGER :: ilatex 1367 CHARACTER(LEN=4) :: Cchar 1368 !#endif 1369 DOUBLE PRECISION :: estimatedTime(1) 1370 1371 IF ( .NOT. this%set ) & 1372 CALL ERROR("Ctqmc_run : QMC not set up ") 1373 IF ( .NOT. this%setU ) & 1374 CALL ERROR("Ctqmc_run : QMC does not have a U this ") 1375 1376 1377 ! OPTIONS of the run 1378 IF ( PRESENT( opt_check ) ) THEN 1379 this%opt_check = opt_check 1380 CALL ImpurityOperator_doCheck(this%Impurity,opt_check) 1381 CALL BathOperator_doCheck(this%Bath,opt_check) 1382 END IF 1383 IF ( PRESENT( opt_movie ) ) & 1384 this%opt_movie = opt_movie 1385 IF ( PRESENT( opt_analysis ) ) & 1386 this%opt_analysis = opt_analysis 1387 IF ( PRESENT ( opt_order ) ) & 1388 this%opt_order = opt_order 1389 IF ( PRESENT ( opt_histo ) ) & 1390 this%opt_histo = opt_histo 1391 IF ( PRESENT ( opt_noise ) ) THEN 1392 this%opt_noise = opt_noise 1393 END IF 1394 IF ( PRESENT ( opt_spectra ) ) & 1395 this%opt_spectra = opt_spectra 1396 1397 this%modGlobalMove(1) = this%sweeps+1 ! No Global Move 1398 this%modGlobalMove(2) = 0 1399 IF ( PRESENT ( opt_gMove ) ) THEN 1400 IF ( opt_gMove .LE. 0 .OR. opt_gMove .GT. this%sweeps ) THEN 1401 this%modGlobalMove(1) = this%sweeps+1 1402 CALL WARNALL("Ctqmc_run : global moves option is <= 0 or > sweeps/cpu -> No global Moves") 1403 ELSE 1404 this%modGlobalMove(1) = opt_gMove 1405 END IF 1406 END IF 1407 1408 CALL Ctqmc_allocateOpt(this) 1409 1410 !#ifdef CTCtqmc_MOVIE 1411 ilatex = 0 1412 IF ( this%opt_movie .EQ. 1 ) THEN 1413 Cchar ="0000" 1414 WRITE(Cchar,'(I4)') this%rank 1415 ilatex = 87+this%rank 1416 OPEN(UNIT=ilatex, FILE="Movie_"//TRIM(ADJUSTL(Cchar))//".tex") 1417 WRITE(ilatex,'(A)') "\documentclass{beamer}" 1418 WRITE(ilatex,'(A)') "\usepackage{color}" 1419 WRITE(ilatex,'(A)') "\setbeamersize{sidebar width left=0pt}" 1420 WRITE(ilatex,'(A)') "\setbeamersize{sidebar width right=0pt}" 1421 WRITE(ilatex,'(A)') "\setbeamersize{text width left=0pt}" 1422 WRITE(ilatex,'(A)') "\setbeamersize{text width right=0pt}" 1423 WRITE(ilatex,*) 1424 WRITE(ilatex,'(A)') "\begin{document}" 1425 WRITE(ilatex,*) 1426 END IF 1427 !#endif 1428 1429 IF ( this%rank .EQ. 0 ) THEN 1430 WRITE(this%ostream,'(A29)') "Starting QMC (Thermalization)" 1431 END IF 1432 1433 !================================= 1434 ! STARTING THERMALIZATION 1435 !================================= 1436 CALL Ctqmc_loop(this,this%thermalization,ilatex) 1437 !================================= 1438 ! ENDING THERMALIZATION 1439 !================================= 1440 1441 estimatedTime = this%runTime 1442 #ifdef HAVE_MPI 1443 CALL MPI_REDUCE([this%runTime], estimatedTime, 1, MPI_DOUBLE_PRECISION, MPI_MAX, & 1444 0, this%MY_COMM, ierr) 1445 #endif 1446 1447 IF ( this%rank .EQ. 0 ) THEN 1448 WRITE(this%ostream,'(A26,I6,A11)') "Thermalization done in ", CEILING(estimatedTime), " seconds" 1449 WRITE(this%ostream,'(A25,I7,A15,I5,A5)') "The QMC should run in ", & 1450 CEILING(estimatedTime(1)*DBLE(this%sweeps)/DBLE(this%thermalization)),& 1451 " seconds on ", this%size, " CPUs" 1452 END IF 1453 1454 !================================= 1455 ! CLEANING CTQMC 1456 !================================= 1457 CALL Ctqmc_clear(this) 1458 1459 !================================= 1460 ! STARTING CTQMC 1461 !================================= 1462 CALL Ctqmc_loop(this,this%sweeps,ilatex) 1463 !================================= 1464 ! ENDING CTQMC 1465 !================================= 1466 1467 IF ( this%opt_movie .EQ. 1 ) THEN 1468 WRITE(ilatex,*) "" 1469 WRITE(ilatex,'(A14)') "\end{document}" 1470 CLOSE(ilatex) 1471 END IF 1472 1473 this%done = .TRUE. 1474 1475 END SUBROUTINE Ctqmc_run
ABINIT/m_Ctqmc/Ctqmc_setG0wTab [ Functions ]
NAME
Ctqmc_setG0wTab
FUNCTION
Set Gow from input array
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc Gomega=G0w opt_fk=F is already inversed with out iwn
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
854 SUBROUTINE Ctqmc_setG0wTab(this,Gomega,opt_fk) 855 856 !Arguments ------------------------------------ 857 TYPE(Ctqmc), INTENT(INOUT) :: this 858 COMPLEX(KIND=8), DIMENSION(:,:), INTENT(IN ) :: Gomega 859 INTEGER , INTENT(IN ) :: opt_fk 860 !Local variable ------------------------------- 861 DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: F 862 863 IF ( .NOT. this%para ) & 864 CALL ERROR("Ctqmc_setG0wTab : Ctqmc_setParameters never called ") 865 866 MALLOC(F,(1:this%samples+1,1:this%flavors)) 867 CALL Ctqmc_computeF(this,Gomega, F, opt_fk) ! mu is changed 868 CALL BathOperator_setF(this%Bath, F) 869 !CALL BathOperator_printF(this%Bath) 870 FREE(F) 871 872 IF ( this%opt_levels .NE. 1 ) THEN ! We compute the mu by hand in computeF 873 CALL ImpurityOperator_setMu(this%Impurity,this%mu) 874 END IF 875 876 this%inF = .TRUE. 877 this%set = .TRUE. 878 879 END SUBROUTINE Ctqmc_setG0wTab
ABINIT/m_Ctqmc/Ctqmc_setMu [ Functions ]
NAME
Ctqmc_setMu
FUNCTION
impose energy levels
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc levels=energy levels vector
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
1115 SUBROUTINE Ctqmc_setMu(this, levels) 1116 1117 !Arguments ------------------------------------ 1118 TYPE(Ctqmc) , INTENT(INOUT) :: this 1119 DOUBLE PRECISION, DIMENSION(:), INTENT(IN ) :: levels 1120 1121 IF ( this%flavors .NE. SIZE(levels,1) ) & 1122 CALL WARNALL("Ctqmc_setMu : Taking energy levels from weiss G(iw)") 1123 1124 this%mu(:)=-levels ! levels = \epsilon_j - \mu 1125 !this%mu =\tilde{\mu} = \mu -\epsilon_j 1126 CALL ImpurityOperator_setMu(this%Impurity,this%mu) 1127 this%opt_levels = 1 1128 END SUBROUTINE Ctqmc_setMu
ABINIT/m_Ctqmc/Ctqmc_setParameters [ Functions ]
NAME
Ctqmc_setParameters
FUNCTION
set all parameters and operators
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc buffer=input parameters
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
475 SUBROUTINE Ctqmc_setParameters(this,buffer) 476 477 !Arguments ------------------------------------ 478 TYPE(Ctqmc), INTENT(INOUT) :: this 479 DOUBLE PRECISION, DIMENSION(1:10), INTENT(IN ) :: buffer 480 481 482 this%thermalization = INT(buffer(3)) !this%thermalization 483 CALL Ctqmc_setSeed(this,INT(buffer(1))) 484 CALL Ctqmc_setSweeps(this,buffer(2)) 485 486 this%measurements = INT(buffer(4)) !this%measurements 487 this%flavors = INT(buffer(5)) !this%flavors 488 this%samples = INT(buffer(6)) !this%samples 489 this%beta = buffer(7) !this%beta 490 this%U = buffer(8) !U 491 this%nspinor = INT(buffer(10))!this%nspinor 492 ! this%mu = buffer(9) !this%mu 493 !this%Wmax = INT(buffer(9)) !Freq 494 !#ifdef CTCtqmc_ANALYSIS 495 ! this%order = INT(buffer(10)) ! order 496 this%inv_dt = this%samples / this%beta 497 !#endif 498 499 !CALL ImpurityOperator_init(this%Impurity,this%flavors,this%beta, this%samples) 500 CALL ImpurityOperator_init(this%Impurity,this%flavors,this%beta) 501 IF ( this%U .GE. 0.d0 ) THEN 502 CALL ImpurityOperator_computeU(this%Impurity,this%U,0.d0) 503 this%setU = .TRUE. 504 END IF 505 ! this%mu = this%mu + this%Impurity%shift_mu 506 507 CALL BathOperator_init(this%Bath, this%flavors, this%samples, this%beta, INT(buffer(9))) 508 509 this%para = .TRUE. 510 511 END SUBROUTINE Ctqmc_setParameters
ABINIT/m_Ctqmc/Ctqmc_setSeed [ Functions ]
NAME
Ctqmc_setSeed
FUNCTION
initialize random number generator
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc iseed=seed from imput
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
594 SUBROUTINE Ctqmc_setSeed(this,iseed) 595 596 !Arguments ------------------------------------ 597 TYPE(Ctqmc), INTENT(INOUT) :: this 598 INTEGER , INTENT(IN ) :: iseed 599 !Local variables ------------------------------ 600 !INTEGER :: n 601 !INTEGER :: i 602 !INTEGER, DIMENSION(:), ALLOCATABLE :: seed 603 604 605 !CALL RANDOM_SEED(size = n) 606 !MALLOC(seed,(n)) 607 !seed = iseed + (/ (i - 1, i = 1, n) /) 608 609 !CALL RANDOM_SEED(PUT = seed+this%rank) 610 611 !FREE(seed) 612 613 this%seed=INT(iseed+this%rank,8) 614 615 END SUBROUTINE Ctqmc_setSeed
ABINIT/m_Ctqmc/Ctqmc_setSweeps [ Functions ]
NAME
Ctqmc_setSweeps
FUNCTION
set the number of sweeps
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc sweeps=asked sweeps
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
539 SUBROUTINE Ctqmc_setSweeps(this,sweeps) 540 541 !Arguments ------------------------------------ 542 TYPE(Ctqmc) , INTENT(INOUT) :: this 543 DOUBLE PRECISION , INTENT(IN ) :: sweeps 544 545 this%sweeps = NINT(sweeps / DBLE(this%size)) 546 ! write(6,*) this%sweeps,NINT(sweeps / DBLE(this%size)),ANINT(sweeps/DBLE(this%size)) 547 IF ( DBLE(this%sweeps) .NE. ANINT(sweeps/DBLE(this%size)) ) & 548 CALL ERROR("Ctqmc_setSweeps : sweeps is negative or too big ") 549 IF ( this%sweeps .LT. 2*CTQMC_SLICE1 ) THEN !202 550 CALL WARNALL("Ctqmc_setSweeps : # sweeps automtically changed ") 551 this%sweeps = 2*CTQMC_SLICE1 552 ! ELSE IF ( this%sweeps .LT. this%thermalization ) THEN 553 ! CALL WARNALL("Ctqmc_setSweeps : Thermalization > sweeps / cpu -> auto fix") 554 ! this%sweeps = this%thermalization 555 END IF 556 IF ( DBLE(NINT(DBLE(this%sweeps)*DBLE(this%size)/DBLE(CTQMC_SLICE1))) .NE. & 557 ANINT(DBLE(this%sweeps)*DBLE(this%size)/DBLE(CTQMC_SLICE1)) ) THEN 558 this%modNoise1 = this%sweeps 559 ELSE 560 this%modNoise1 = MIN(this%sweeps,INT(DBLE(this%sweeps)*DBLE(this%size) / DBLE(CTQMC_SLICE1))) !101 561 END IF 562 this%modNoise2 = MAX(this%modNoise1 / CTQMC_SLICE2, 1) ! 100 563 ! this%modGlobalMove(1) = this%thermalization / 10 + 1 564 ! this%modGlobalMove(2) = 0 565 566 END SUBROUTINE Ctqmc_setSweeps
ABINIT/m_Ctqmc/Ctqmc_setU [ Functions ]
NAME
Ctqmc_setU
FUNCTION
set the interaction this
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc matU=interaction this
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
938 SUBROUTINE Ctqmc_setU(this,matU) 939 940 !Arguments ------------------------------------ 941 TYPE(Ctqmc), INTENT(INOUT) :: this 942 !Local variables ------------------------------ 943 DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: matU 944 945 IF ( SIZE(matU) .NE. this%flavors*this%flavors ) & 946 CALL ERROR("Ctqmc_setU : Wrong interaction this (size) ") 947 948 CALL ImpurityOperator_setUmat(this%Impurity, matU) 949 this%setU = .TRUE. 950 END SUBROUTINE Ctqmc_setU
ABINIT/m_Ctqmc/Ctqmc_symmetrizeGreen [ Functions ]
NAME
Ctqmc_symmetrizeGreen
FUNCTION
optionnaly symmetrize the green functions
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc syms=weight factors
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
2792 SUBROUTINE Ctqmc_symmetrizeGreen(this, syms) 2793 2794 !Arguments ------------------------------------ 2795 TYPE(Ctqmc) , INTENT(INOUT) :: this 2796 DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN ) :: syms 2797 !Local variables ------------------------------ 2798 INTEGER :: iflavor1 2799 INTEGER :: iflavor2 2800 INTEGER :: flavors 2801 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: green_tmp 2802 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(: ) :: n_tmp 2803 2804 flavors = this%flavors 2805 IF ( SIZE(syms,1) .NE. flavors .OR. SIZE(syms,2) .NE. flavors ) THEN 2806 CALL WARNALL("Ctqmc_symmetrizeGreen : wrong opt_sym -> not symmetrizing") 2807 RETURN 2808 END IF 2809 2810 MALLOC(green_tmp,(1:this%samples+1,flavors)) 2811 green_tmp(:,:) = 0.d0 2812 MALLOC(n_tmp,(1:flavors)) 2813 n_tmp(:) = 0.d0 2814 DO iflavor1=1, flavors 2815 DO iflavor2=1,flavors 2816 green_tmp(:,iflavor1) = green_tmp(:,iflavor1) & 2817 + syms(iflavor2,iflavor1) * this%Greens(iflavor2)%oper(:) 2818 n_tmp(iflavor1) = n_tmp(iflavor1) & 2819 + syms(iflavor2,iflavor1) * this%measN(1,iflavor2) 2820 END DO 2821 END DO 2822 DO iflavor1=1, flavors 2823 this%Greens(iflavor1)%oper(:) = green_tmp(:,iflavor1) 2824 this%measN(1,iflavor1) = n_tmp(iflavor1) 2825 END DO 2826 FREE(green_tmp) 2827 FREE(n_tmp) 2828 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
updated=something changed
SIDE EFFECTS
NOTES
SOURCE
1746 SUBROUTINE Ctqmc_tryAddRemove(this,updated) 1747 1748 !Arguments ------------------------------------ 1749 TYPE(Ctqmc) , INTENT(INOUT) :: this 1750 ! TYPE(BathOperator) , INTENT(INOUT) :: Bath 1751 ! TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity 1752 LOGICAL , INTENT( OUT) :: updated 1753 !Local variables ------------------------------ 1754 INTEGER :: position 1755 INTEGER , DIMENSION(1:2) :: nature ! -2 for antiseg and 1 for seg 1756 INTEGER :: i! -2 for antiseg and 1 for seg 1757 DOUBLE PRECISION :: action 1758 DOUBLE PRECISION :: beta 1759 DOUBLE PRECISION :: time1 1760 DOUBLE PRECISION :: time2 1761 DOUBLE PRECISION :: time_avail 1762 DOUBLE PRECISION :: det_ratio 1763 DOUBLE PRECISION :: imp_trace 1764 DOUBLE PRECISION :: signe 1765 DOUBLE PRECISION :: tail 1766 DOUBLE PRECISION, DIMENSION(1:2) :: CdagC_1 1767 1768 IF ( .NOT. this%set ) & 1769 CALL ERROR("Ctqmc_trySegment : QMC not set ") 1770 1771 nature(1) = CTQMC_SEGME 1772 nature(2) = CTQMC_ANTIS 1773 beta = this%beta 1774 1775 updated = .FALSE. 1776 tail = DBLE(this%Impurity%particles(this%Impurity%activeFlavor)%tail) 1777 1778 1779 DO i = 1, 2 1780 signe = SIGN(1.d0,DBLE(nature(i))) 1781 1782 !CALL RANDOM_NUMBER(action) 1783 CALL OurRng(this%seed,action) 1784 1785 IF ( action .LT. .5d0 ) THEN ! Ajout de segment 1786 !CALL RANDOM_NUMBER(time1) 1787 CALL OurRng(this%seed,time1) 1788 time1 = time1 * beta 1789 time_avail = ImpurityOperator_getAvailableTime(this%Impurity,time1,position) * signe 1790 IF ( time_avail .GT. 0.d0 ) THEN 1791 !CALL RANDOM_NUMBER(time2) 1792 CALL OurRng(this%seed,time2) 1793 IF ( time2 .EQ. 0.d0 ) THEN 1794 CALL OurRng(this%seed,time2) ! Prevent null segment 1795 END IF 1796 time2 = time1 + time2 * time_avail 1797 CdagC_1(Cdag_) = ((1.d0+signe)*time1+(1.d0-signe)*time2)*0.5d0 1798 CdagC_1(C_ ) = ((1.d0+signe)*time2+(1.d0-signe)*time1)*0.5d0 1799 det_ratio = BathOperator_getDetAdd(this%Bath,CdagC_1,position,this%Impurity%particles(this%Impurity%activeFlavor)) 1800 imp_trace = ImpurityOperator_getTraceAdd(this%Impurity,CdagC_1) 1801 !CALL RANDOM_NUMBER(time1) 1802 CALL OurRng(this%seed,time1) 1803 IF ( det_ratio*imp_trace .LT. 0.d0 ) THEN 1804 this%stats(nature(i)+CTQMC_DETSI) = this%stats(nature(i)+CTQMC_DETSI) + 1.d0 1805 END IF 1806 IF ( (time1 * (tail + 1.d0 )) & 1807 .LT. (beta * time_avail * det_ratio * imp_trace ) ) THEN 1808 CALL ImpurityOperator_add(this%Impurity,CdagC_1,position) 1809 CALL BathOperator_setMAdd(this%bath,this%Impurity%particles(this%Impurity%activeFlavor)) 1810 this%stats(nature(i)+CTQMC_ADDED) = this%stats(nature(i)+CTQMC_ADDED) + 1.d0 1811 updated = .TRUE. .OR. updated 1812 tail = tail + 1.d0 1813 END IF 1814 END IF 1815 1816 ELSE ! Supprimer un segment 1817 IF ( tail .GT. 0.d0 ) THEN 1818 !CALL RANDOM_NUMBER(time1) 1819 CALL OurRng(this%seed,time1) 1820 position = INT(((time1 * tail) + 1.d0) * signe ) 1821 time_avail = ImpurityOperator_getAvailedTime(this%Impurity,position) 1822 det_ratio = BathOperator_getDetRemove(this%Bath,position) 1823 imp_trace = ImpurityOperator_getTraceRemove(this%Impurity,position) 1824 !CALL RANDOM_NUMBER(time1) 1825 CALL OurRng(this%seed,time1) 1826 IF ( det_ratio * imp_trace .LT. 0.d0 ) THEN 1827 this%stats(nature(i)+CTQMC_DETSI) = this%stats(nature(i)+CTQMC_DETSI) + 1.d0 1828 END IF 1829 IF ( (time1 * beta * time_avail ) & 1830 .LT. (tail * det_ratio * imp_trace) ) THEN 1831 CALL ImpurityOperator_remove(this%Impurity,position) 1832 CALL BathOperator_setMRemove(this%Bath,this%Impurity%particles(this%Impurity%activeFlavor)) 1833 this%stats(nature(i)+CTQMC_REMOV) = this%stats(nature(i)+CTQMC_REMOV) + 1.d0 1834 updated = .TRUE. .OR. updated 1835 tail = tail -1.d0 1836 END IF 1837 END IF 1838 END IF 1839 END DO 1840 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-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
this=ctqmc
OUTPUT
flav_i=first flavor swaped flav_j=second flavor swaped
SIDE EFFECTS
NOTES
SOURCE
1869 SUBROUTINE Ctqmc_trySwap(this,flav_i,flav_j) 1870 1871 !Arguments ------------------------------------ 1872 TYPE(Ctqmc) , INTENT(INOUT) :: this 1873 ! TYPE(BathOperator) , INTENT(INOUT) :: Bath 1874 ! TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity 1875 INTEGER , INTENT( OUT) :: flav_i 1876 INTEGER , INTENT( OUT) :: flav_j 1877 !Local variables ------------------------------ 1878 INTEGER :: flavor_i 1879 INTEGER :: flavor_j 1880 DOUBLE PRECISION :: rnd 1881 DOUBLE PRECISION :: lengthi 1882 DOUBLE PRECISION :: lengthj 1883 DOUBLE PRECISION :: overlapic1 1884 DOUBLE PRECISION :: overlapjc1 1885 DOUBLE PRECISION :: overlapic2 1886 DOUBLE PRECISION :: overlapjc2 1887 DOUBLE PRECISION :: detic1 1888 DOUBLE PRECISION :: detjc1 1889 DOUBLE PRECISION :: detic2 1890 DOUBLE PRECISION :: detjc2 1891 DOUBLE PRECISION :: det_ratio 1892 DOUBLE PRECISION :: local_ratio 1893 1894 !CALL RANDOM_NUMBER(rnd) 1895 CALL OurRng(this%seed,rnd) 1896 flavor_i = NINT(rnd*DBLE(this%flavors-1.d0))+1 1897 !CALL RANDOM_NUMBER(rnd) 1898 CALL OurRng(this%seed,rnd) 1899 flavor_j = NINT(rnd*DBLE(this%flavors-1.d0))+1 1900 1901 flav_i = 0 1902 flav_j = 0 1903 1904 IF ( flavor_i .NE. flavor_j ) THEN 1905 ! On tente d'intervertir i et j 1906 ! Configuration actuelle : 1907 this%modGlobalMove(2) = this%modGlobalMove(2)+1 1908 detic1 = BathOperator_getDetF(this%Bath,flavor_i) 1909 detjc1 = BathOperator_getDetF(this%Bath,flavor_j) 1910 lengthi = ImpurityOperator_measN(this%Impurity,flavor_i) 1911 lengthj = ImpurityOperator_measN(this%Impurity,flavor_j) 1912 overlapic1 = ImpurityOperator_overlapFlavor(this%Impurity,flavor_i) 1913 overlapjc1 = ImpurityOperator_overlapFlavor(this%Impurity,flavor_j) 1914 ! Configuration nouvelle : 1915 detic2 = BathOperator_getDetF(this%Bath,flavor_i,this%Impurity%particles(flavor_j)) 1916 detjc2 = BathOperator_getDetF(this%Bath,flavor_j,this%Impurity%particles(flavor_i)) 1917 ! lengths unchanged 1918 overlapic2 = ImpurityOperator_overlapSwap(this%Impurity,flavor_i,flavor_j) 1919 overlapjc2 = ImpurityOperator_overlapSwap(this%Impurity,flavor_j,flavor_i) 1920 1921 ! IF ( detic1*detjc1 .EQ. detic2*detjc2 ) THEN 1922 ! det_ratio = 1.d0 1923 ! ELSE IF ( detic1*detjc1 .EQ. 0.d0 ) THEN 1924 ! det_ratio = detic2*detjc2 ! evite de diviser par 0 si pas de segment 1925 ! ELSE 1926 det_ratio = detic2*detjc2/(detic1*detjc1) 1927 ! END IF 1928 local_ratio = DEXP(-overlapic2*overlapjc2+overlapic1*overlapjc1 & 1929 +(lengthj-lengthi)*(this%mu(flavor_i)-this%mu(flavor_j))) 1930 ! Wloc = exp(muN-Uo) 1931 !CALL RANDOM_NUMBER(rnd) 1932 CALL OurRng(this%seed,rnd) 1933 IF ( rnd .LT. local_ratio*det_ratio ) THEN ! swap accepted 1934 CALL ImpurityOperator_swap(this%Impurity, flavor_i,flavor_j) 1935 CALL BathOperator_swap (this%Bath , flavor_i,flavor_j) 1936 this%swap = this%swap + 1.d0 1937 flav_i = flavor_i 1938 flav_j = flavor_j 1939 ! ELSE 1940 ! CALL WARN("Swap refused") 1941 ! WRITE(this%ostream,'(6E24.14)') local_ratio, det_ratio, detic1, detjc1, detic2, detjc2 1942 END IF 1943 END IF 1944 1945 END SUBROUTINE Ctqmc_trySwap
ABINIT/m_Ctqmcoffdiag/Ctqmc_setMagmom [ Functions ]
NAME
Ctqmcoffdiag_setMagmom
FUNCTION
set the Magnetic moment matrix for susceptibility
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (F. Gendron) 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
OUTPUT
SIDE EFFECTS
NOTES
PARENTS
Will be filled automatically by the parent script
CHILDREN
Will be filled automatically by the parent script
SOURCE
3683 SUBROUTINE Ctqmc_setMagmom(this,Magmom_orb,Magmom_spin,Magmom_tot) 3684 3685 !Arguments ------------------------------------ 3686 TYPE(Ctqmc), INTENT(INOUT) :: this 3687 DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: Magmom_orb 3688 DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: Magmom_spin 3689 DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: Magmom_tot 3690 !Local variables ------------------------------ 3691 ! INTEGER :: iflavor1,iflavor2 3692 3693 ! do iflavor1=1,10 3694 ! do iflavor2=1,10 3695 ! if(iflavor1==iflavor2) THEN 3696 ! write(6,*) iflavor1, iflavor2, Magmom(iflavor1,iflavor2) 3697 ! end if 3698 ! end do 3699 ! end do 3700 3701 IF ( SIZE(Magmom_orb) .NE. this%flavors*this%flavors ) & 3702 CALL ERROR("Ctqmc_setMagmomm : Wrong Magnetic Moment matrix (size) ") 3703 3704 CALL ImpurityOperator_setMagmommat(this%Impurity, Magmom_orb, Magmom_spin, Magmom_tot) 3705 3706 END SUBROUTINE Ctqmc_setMagmom
m_Ctqmc/Ctqmc [ Types ]
NAME
Ctqmc
FUNCTION
This structured datatype contains the necessary data
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
74 TYPE, PUBLIC :: Ctqmc 75 76 LOGICAL _PRIVATE :: init = .FALSE. 77 ! Flag: is MC initialized 78 79 LOGICAL _PRIVATE :: set = .FALSE. 80 ! Flag: ?? 81 82 LOGICAL _PRIVATE :: setU = .FALSE. 83 ! Flag: is U Set ? 84 85 LOGICAL _PRIVATE :: inF = .FALSE. 86 ! Flag: is hybridization fct in input ? 87 88 LOGICAL _PRIVATE :: done = .FALSE. 89 ! Flag: is MC terminated ? 90 91 LOGICAL _PRIVATE :: para = .FALSE. 92 ! Flag: do we have parameters in input 93 94 LOGICAL _PRIVATE :: have_MPI = .FALSE. 95 ! Flag: 96 97 INTEGER _PRIVATE :: opt_movie = 0 98 ! 99 100 INTEGER _PRIVATE :: opt_analysis = 0 101 ! correlations 102 103 INTEGER _PRIVATE :: opt_check = 0 104 ! various check 0 105 ! various check 1 impurity 106 ! various check 2 bath 107 ! various check 3 both 108 109 INTEGER _PRIVATE :: opt_order = 0 110 ! nb of segments max for analysis 111 112 INTEGER _PRIVATE :: opt_histo = 0 113 ! Enable histograms 114 115 INTEGER _PRIVATE :: opt_noise = 0 116 ! compute noise 117 118 INTEGER _PRIVATE :: opt_spectra = 0 119 ! markov chain FT (correlation time) 120 121 INTEGER _PRIVATE :: opt_levels = 0 122 ! do we have energy levels 123 124 INTEGER _PRIVATE :: flavors 125 ! 126 INTEGER _PRIVATE :: nspinor 127 ! 128 INTEGER _PRIVATE :: measurements 129 ! nb of measure in the MC 130 131 INTEGER _PRIVATE :: samples 132 ! nb of L points 133 134 INTEGER(8) _PRIVATE :: seed 135 ! 136 137 INTEGER _PRIVATE :: sweeps 138 ! 139 140 INTEGER _PRIVATE :: thermalization 141 ! 142 143 INTEGER _PRIVATE :: ostream 144 ! output file 145 146 INTEGER _PRIVATE :: istream 147 ! input file 148 149 INTEGER _PRIVATE :: modNoise1 150 ! measure the noise each modNoise1 151 152 INTEGER _PRIVATE :: modNoise2 153 ! measure the noise each modNoise2 154 155 INTEGER _PRIVATE :: activeFlavor 156 ! orbital on which one do sth now 157 158 INTEGER, DIMENSION(1:2) _PRIVATE :: modGlobalMove 159 ! 1: gloabl move each modglobalmove(1) 160 ! 2: we have done modglobalmove(2) for two different orbitals. 161 162 INTEGER _PRIVATE :: Wmax 163 ! Max freq for FT 164 165 DOUBLE PRECISION, DIMENSION(1:6) _PRIVATE :: stats 166 ! to now how many negative determinant, antisegments,seeme.e.twfs...j 167 168 DOUBLE PRECISION _PRIVATE :: swap 169 ! nb of successfull GM 170 171 INTEGER _PRIVATE :: MY_COMM 172 ! 173 174 INTEGER _PRIVATE :: rank 175 ! 176 177 INTEGER _PRIVATE :: size 178 ! size of MY_COMM 179 180 DOUBLE PRECISION _PRIVATE :: runTime ! time for the run routine 181 ! 182 183 DOUBLE PRECISION _PRIVATE :: beta 184 ! 185 186 DOUBLE PRECISION _PRIVATE :: U 187 188 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) _PRIVATE :: mu 189 ! levels 190 191 TYPE(GreenHyb) , ALLOCATABLE, DIMENSION(: ) _PRIVATE :: Greens 192 193 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) _PRIVATE :: measN 194 ! measure of occupations (3or4,flavor) 195 196 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) _PRIVATE :: measDE 197 ! (flavor,flavor) double occupancies 198 ! (1,1): total energy of correlation. 199 200 DOUBLE PRECISION _PRIVATE :: a_Noise 201 ! Noise a exp (-bx) for the noise 202 203 DOUBLE PRECISION _PRIVATE :: b_Noise 204 ! Noise a exp (-bx) for the noise 205 206 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: abNoiseG !(ab,tau,flavor) 207 ! Noise but for G 208 209 TYPE(Vector) , DIMENSION(1:2) _PRIVATE :: measNoise 210 TYPE(Vector), ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: measNoiseG !(tau,flavor,mod) 211 ! accumulate each value relataed to measurenoise 1 2 212 213 DOUBLE PRECISION _PRIVATE :: inv_dt 214 ! 1/(beta/L) 215 216 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,: ) _PRIVATE :: measPerturbation 217 ! opt_order,nflavor 218 219 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: occup_histo_time 220 ! nflavor 221 222 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: occupconfig 223 ! 2**nflavor 224 225 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: suscep 226 ! samples 227 228 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: chi 229 ! samples 230 231 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: chicharge 232 ! samples 233 234 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: ntot 235 ! occupation total, t2g, eg 236 237 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: measCorrelation 238 ! segment,antisegment,nflavor,nflavor 239 240 DOUBLE PRECISION _PRIVATE :: errorImpurity 241 ! check 242 243 DOUBLE PRECISION _PRIVATE :: errorBath 244 ! for check 245 246 TYPE(BathOperator) _PRIVATE :: Bath 247 248 TYPE(ImpurityOperator) _PRIVATE :: Impurity 249 250 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) _PRIVATE :: density 251 252 END TYPE Ctqmc