TABLE OF CONTENTS


ABINIT/m_Ctqmc [ Modules ]

[ Top ] [ Modules ]

NAME

  m_Ctqmc

FUNCTION

  Manage and drive all the CTQMC
  Should not be used if you don't know what you do
  Please use CtqmcInterface

COPYRIGHT

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

NOTES

SOURCE

25 #include "defs.h"
26 
27 MODULE m_Ctqmc
28 
29 USE m_Global
30 USE m_GreenHyb
31 USE m_BathOperator
32 USE m_ImpurityOperator
33 USE m_Stat
34 USE m_FFTHyb
35 USE m_OurRng
36 USE m_Vector
37 #ifdef HAVE_MPI2
38 USE mpi
39 #endif
40 
41 IMPLICIT NONE

ABINIT/m_Ctqmc/Ctqmc_allocateAll [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_allocateAll

FUNCTION

  Allocate all non option varibales

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

622 SUBROUTINE Ctqmc_allocateAll(this)
623 
624 !Arguments ------------------------------------
625   TYPE(Ctqmc), INTENT(INOUT) :: this
626 !Local variables ------------------------------
627   INTEGER                  :: flavors
628 
629   IF ( .NOT. this%para ) &
630     CALL ERROR("Ctqmc_allocateAll : Ctqmc_setParameters never called  ")
631 
632   flavors = this%flavors
633 
634   DT_FREEIF(this%Greens)
635   DT_MALLOC(this%Greens,(1:flavors))
636 
637   FREEIF(this%measN)
638   MALLOC(this%measN,(1:4,1:flavors))
639   this%measN = 0.d0
640 
641   FREEIF(this%measDE)
642   MALLOC(this%measDE,(1:flavors,1:flavors) )
643   this%measDE = 0.d0
644 
645   FREEIF(this%mu)
646   MALLOC(this%mu,(1:flavors) )
647   this%mu = 0.d0
648 END SUBROUTINE Ctqmc_allocateAll

ABINIT/m_Ctqmc/Ctqmc_allocateOpt [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_allocateOpt

FUNCTION

  allocate all option variables 

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

676 SUBROUTINE Ctqmc_allocateOpt(this)
677 
678 !Arguments ------------------------------------
679   TYPE(Ctqmc), INTENT(INOUT) :: this
680 !Local variables ------------------------------
681   INTEGER :: i
682   INTEGER :: j
683   INTEGER :: k
684 
685   IF ( .NOT. this%para ) &
686     CALL ERROR("Ctqmc_allocateOpt : Ctqmc_setParameters never called  ")
687 
688   IF ( this%opt_analysis .EQ. 1 ) THEN
689     FREEIF(this%measCorrelation)
690     MALLOC(this%measCorrelation,(1:this%samples+1,1:3,1:this%flavors))
691     this%measCorrelation = 0.d0
692   END IF
693 
694   IF ( this%opt_order .GT. 0 ) THEN
695     FREEIF(this%measPerturbation)
696     MALLOC(this%measPerturbation,(1:this%opt_order,1:this%flavors))
697     this%measPerturbation = 0.d0
698   END IF
699 
700   IF ( this%opt_histo .GT. 0 ) THEN
701     FREEIF(this%occup_histo_time)
702     MALLOC(this%occup_histo_time,(1:this%flavors+1))
703     this%occup_histo_time= 0.d0
704   END IF
705 
706   IF ( this%opt_noise .EQ. 1 ) THEN
707     IF ( ALLOCATED(this%measNoiseG) ) THEN
708       DO i=1,2
709         DO j = 1, this%flavors
710           DO k= 1, this%samples+1
711             CALL Vector_destroy(this%measNoiseG(k,j,i))
712           END DO
713         END DO
714       END DO
715       DT_FREE(this%measNoiseG)
716     END IF
717     DT_MALLOC(this%measNoiseG,(1:this%samples+1,1:this%flavors,1:2))
718     !DO i=1,2
719       DO j = 1, this%flavors
720         DO k= 1, this%samples+1
721           CALL Vector_init(this%measNoiseG(k,j,1),CTQMC_SLICE1)
722         END DO
723       END DO
724       DO j = 1, this%flavors
725         DO k= 1, this%samples+1
726           CALL Vector_init(this%measNoiseG(k,j,2),CTQMC_SLICE1*CTQMC_SLICE2+1) ! +1 pour etre remplacer ceil
727         END DO
728       END DO
729     !END DO
730     FREEIF(this%abNoiseG)
731     MALLOC(this%aBNoiseG,(1:2,1:this%samples+1,this%flavors))
732     this%abNoiseG = 0.d0
733   END IF
734 
735   IF (this%opt_spectra .GE. 1 ) THEN
736     FREEIF(this%density)
737     !MALLOC(this%density,(1:this%thermalization,1:this%flavors))
738     i = CEILING(DBLE(this%thermalization+this%sweeps)/DBLE(this%measurements*this%opt_spectra))
739     MALLOC(this%density,(1:this%flavors+1,1:i))
740     this%density = 0.d0
741   END IF
742 !#endif
743 END SUBROUTINE Ctqmc_allocateOpt

ABINIT/m_Ctqmc/Ctqmc_clear [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_clear

FUNCTION

  clear a ctqmc run

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

937 SUBROUTINE Ctqmc_clear(this)
938 
939 !Arguments ------------------------------------
940   TYPE(Ctqmc), INTENT(INOUT) :: this
941 !Local variables ------------------------------
942   INTEGER :: i
943   INTEGER :: j
944   INTEGER :: k
945 
946   this%measN(1,:) = 0.d0
947   this%measN(2,:) = 0.d0
948   !Do not set measN(3,:) to 0 to avoid erasing N between therm and ctqmc
949   this%measN(4,:) = 0.d0
950   this%measDE = 0.d0
951 !  this%seg_added    = 0.d0
952 !  this%anti_added   = 0.d0
953 !  this%seg_removed  = 0.d0
954 !  this%anti_removed = 0.d0
955 !  this%seg_sign     = 0.d0
956 !  this%anti_sign    = 0.d0
957   this%stats(:)     = 0.d0
958   this%swap         = 0.d0
959   this%runTime      = 0.d0
960   this%modGlobalMove(2) = 0 
961   CALL Vector_clear(this%measNoise(1))
962   CALL Vector_clear(this%measNoise(2))
963 !#ifdef CTCtqmc_CHECK
964   this%errorImpurity = 0.d0
965   this%errorBath     = 0.d0
966 !#endif
967   DO j = 1, this%flavors
968     CALL GreenHyb_clear(this%Greens(j))
969   END DO
970 !#ifdef CTCtqmc_ANALYSIS
971   IF ( this%opt_analysis .EQ. 1 .AND. ALLOCATED(this%measCorrelation) ) &    
972     this%measCorrelation = 0.d0 
973   IF ( this%opt_order .GT. 0 .AND. ALLOCATED(this%measPerturbation) ) &
974     this%measPerturbation = 0.d0
975   IF ( this%opt_noise .EQ. 1 .AND. ALLOCATED(this%measNoiseG) ) THEN
976     DO i=1,2
977       DO j = 1, this%flavors
978         DO k= 1, this%samples+1
979           CALL Vector_clear(this%measNoiseG(k,j,i))
980         END DO
981       END DO
982     END DO
983   END IF
984 !#endif
985 END SUBROUTINE Ctqmc_clear

ABINIT/m_Ctqmc/Ctqmc_computeF [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_computeF

FUNCTION

  Compute the hybridization function

COPYRIGHT

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

INPUTS

  this=ctqmc
  Gomega=G0 to compute F
  opt_fk=What is Gomega

OUTPUT

  F=hybridization function

SIDE EFFECTS

NOTES

SOURCE

1118 SUBROUTINE Ctqmc_computeF(this, Gomega, F, opt_fk)
1119 
1120 !Arguments ------------------------------------
1121   TYPE(Ctqmc)                       , INTENT(INOUT) :: this
1122   COMPLEX(KIND=8), DIMENSION(:,:), INTENT(IN   ) :: Gomega
1123   !INTEGER                         , INTENT(IN   ) :: Wmax
1124   DOUBLE PRECISION, DIMENSION(:,:), INTENT(INOUT) :: F
1125   INTEGER                         , INTENT(IN   ) :: opt_fk
1126 !Local variables ------------------------------
1127   INTEGER                                         :: flavors
1128   INTEGER                                         :: samples
1129   INTEGER                                         :: iflavor
1130   INTEGER                                         :: iomega
1131   INTEGER                                         :: itau
1132   DOUBLE PRECISION                                :: pi_invBeta
1133   DOUBLE PRECISION                                :: K
1134   DOUBLE PRECISION                                :: re
1135   DOUBLE PRECISION                                :: im
1136   COMPLEX(KIND=8), DIMENSION(:,:), ALLOCATABLE   :: F_omega
1137   TYPE(GreenHyb)                                     :: F_tmp
1138 
1139   flavors    = this%flavors
1140 
1141   samples    = this%samples
1142   pi_invBeta = ACOS(-1.d0) / this%beta
1143   this%Wmax=SIZE(Gomega,1)
1144 
1145   IF ( this%have_MPI .EQV. .TRUE. ) THEN
1146     CALL GreenHyb_init(F_tmp,samples,this%beta,MY_COMM=this%MY_COMM)
1147   ELSE
1148     CALL GreenHyb_init(F_tmp,samples,this%beta)
1149   END IF
1150 !  K = this%mu
1151 
1152   MALLOC(F_omega,(1:this%Wmax,1:flavors))
1153 
1154   !IF ( this%rank .EQ. 0 ) &
1155     !OPEN(UNIT=9876,FILE="K.dat",POSITION="APPEND")
1156   IF ( opt_fk .EQ. 0 ) THEN
1157     DO iflavor = 1, flavors
1158       DO iomega=1,this%Wmax
1159         re = REAL(Gomega(iomega,iflavor))
1160         im = AIMAG(Gomega(iomega,iflavor))
1161         F_omega(iomega,iflavor) = CMPLX(-re/(re*re+im*im),im/(re*re+im*im),8)
1162       END DO
1163     END DO
1164     !F_omega = CMPLX(-1.d0,0,8)/Gomega
1165   ELSE
1166     F_omega = Gomega
1167   END IF
1168 
1169   DO iflavor = 1, flavors
1170     IF ( this%opt_levels .EQ. 1 ) THEN
1171       K = this%mu(iflavor)
1172     ELSE
1173       K = -REAL(F_omega(this%Wmax, iflavor))
1174 !    this%mu = K
1175       this%mu(iflavor) = K 
1176     END IF
1177     !IF ( this%rank .EQ. 0 ) &
1178     !WRITE(9876,'(I4,2E22.14)') iflavor, K, REAL(-F_omega(this%Wmax, iflavor))
1179     !IF(this%rank .EQ.0) &
1180     !WRITE(this%ostream,*) "CTQMC K, this%mu = ",K,this%mu(iflavor)
1181     !WRITE(this%ostream,*) "CTQMC beta     = ",this%beta
1182     IF ( opt_fk .EQ. 0 ) THEN
1183       DO iomega = 1, this%Wmax
1184         re = REAL(F_omega(iomega,iflavor))
1185         im = AIMAG(F_omega(iomega,iflavor))
1186         F_omega(iomega,iflavor) = CMPLX(re + K, im + (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, 8)
1187         !if(iflavor==1.and.this%rank==0) then
1188           !write(224,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(F_omega(iomega,iflavor)),imag(F_omega(iomega,iflavor))
1189           !write(225,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(Gomega(iomega, iflavor)),imag(Gomega(iomega, iflavor))
1190         !end if 
1191       END DO
1192     ELSE
1193       DO iomega = 1, this%Wmax
1194         F_omega(iomega,iflavor) = F_omega(iomega,iflavor) &
1195                     + CMPLX(K, 0.d0, 8)
1196         !if(iflavor==1.and.this%rank==0) then
1197           !write(224,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(F_omega(iomega,iflavor)),imag(F_omega(iomega,iflavor))
1198           !write(225,*) (2.d0*DBLE(iomega)-1.d0) * pi_invBeta, real(Gomega(iomega, iflavor)),imag(Gomega(iomega, iflavor))
1199         !end if 
1200       END DO
1201     END IF
1202     K = REAL(CMPLX(0,(2.d0*DBLE(this%Wmax)-1.d0)*pi_invBeta,8)*F_omega(this%Wmax,iflavor))
1203     CALL GreenHyb_setMuD1(this%Greens(iflavor),this%mu(iflavor),K)
1204     CALL GreenHyb_setOperW(F_tmp,F_omega(:,iflavor))
1205     !CALL GreenHyb_backFourier(F_tmp,F_omega(:,iflavor))
1206     CALL GreenHyb_backFourier(F_tmp)
1207     F(1:samples+1,iflavor) = (/ (-F_tmp%oper(samples+1-itau),itau=0,samples) /)
1208   END DO
1209   IF ( this%rank .EQ. 0 ) THEN
1210     open (unit=346,file='Hybridization.dat',status='unknown',form='formatted')
1211     DO iflavor = 1, flavors
1212       write(346,*) "#",iflavor
1213       do  itau=1,this%samples+1
1214         write(346,*) itau,F(itau,iflavor)
1215       enddo
1216       write(346,*) 
1217     END DO
1218     close(346)
1219   ENDIF
1220   FREE(F_omega)
1221   CALL GreenHyb_destroy(F_tmp)
1222 END SUBROUTINE Ctqmc_computeF

ABINIT/m_Ctqmc/Ctqmc_destroy [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_destroy

FUNCTION

  destroy and deallocate all variables

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

3349 SUBROUTINE Ctqmc_destroy(this)
3350 
3351 !Arguments ------------------------------------
3352   TYPE(Ctqmc), INTENT(INOUT) :: this
3353 !Local variables ------------------------------
3354   INTEGER                  :: iflavor
3355   INTEGER                  :: flavors
3356   INTEGER                  :: i
3357   INTEGER                  :: j
3358   INTEGER                  :: k
3359 
3360   if ( this%init .EQV. .FALSE. ) RETURN
3361 
3362   flavors = this%flavors
3363 
3364   CALL ImpurityOperator_destroy(this%Impurity)
3365   CALL BathOperator_destroy(this%Bath)
3366   CALL Vector_destroy(this%measNoise(1))
3367   CALL Vector_destroy(this%measNoise(2))
3368 
3369   IF ( ALLOCATED(this%Greens) ) THEN
3370     DO iflavor = 1, flavors
3371      CALL GreenHyb_destroy(this%Greens(iflavor))
3372     END DO
3373     DT_FREE( this%Greens )
3374   END IF
3375 !#ifdef CTCtqmc_ANALYSIS
3376   FREEIF(this%measCorrelation)
3377   FREEIF(this%measPerturbation)
3378   IF ( this%opt_histo .GT. 0 ) THEN
3379     FREEIF(this%occup_histo_time)
3380   ENDIF
3381   FREEIF(this%measN)
3382   FREEIF(this%measDE)
3383   FREEIF(this%mu)
3384   FREEIF(this%abNoiseG)
3385   IF ( ALLOCATED(this%measNoiseG) ) THEN
3386     DO i=1,2
3387       DO j = 1, this%flavors
3388         DO k= 1, this%samples+1
3389           CALL Vector_destroy(this%measNoiseG(k,j,i))
3390         END DO
3391       END DO
3392     END DO
3393     DT_FREE(this%measNoiseG)
3394   END IF
3395   FREEIF(this%density)
3396 !#endif
3397   this%ostream        = 0
3398   this%istream        = 0
3399  
3400   this%sweeps         = 0
3401   this%thermalization = 0
3402   this%flavors        = 0
3403   this%samples        = 0
3404   this%beta           = 0.d0
3405 !  this%seg_added      = 0.d0
3406 !  this%anti_added     = 0.d0
3407 !  this%seg_removed    = 0.d0
3408 !  this%anti_removed   = 0.d0
3409 !  this%seg_sign       = 0.d0
3410 !  this%anti_sign      = 0.d0
3411   this%stats          = 0.d0
3412   this%swap           = 0.d0
3413 
3414 
3415   this%set  = .FALSE.
3416   this%done = .FALSE.
3417   this%init = .FALSE.
3418 END SUBROUTINE Ctqmc_destroy

ABINIT/m_Ctqmc/Ctqmc_getD [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_getD

FUNCTION

  get double occupation

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

  D=full double occupation

SIDE EFFECTS

NOTES

SOURCE

2710 SUBROUTINE Ctqmc_getD(this, D)
2711 
2712 !Arguments ------------------------------------
2713   TYPE(Ctqmc)       , INTENT(IN ) :: this
2714   DOUBLE PRECISION, DIMENSION(:,:), INTENT(OUT) :: D
2715 !Local variables ------------------------------
2716   INTEGER                       :: iflavor1
2717   INTEGER                       :: iflavor2
2718 
2719   IF ( SIZE(D,1) .LT. this%flavors .OR. SIZE(D,2) .LT. this%flavors ) &
2720     CALL ERROR("Ctqmc_getD : Dimensions of array D are too small")
2721 
2722   D = 0.d0
2723 
2724   DO iflavor1 = 1, this%flavors
2725     DO iflavor2 = 1, this%flavors
2726       D(iflavor2,iflavor1) =  this%measDE(iflavor2,iflavor1)
2727     END DO
2728     D(iflavor1,iflavor1) = 0.d0
2729   END DO
2730 
2731 END SUBROUTINE Ctqmc_getD

ABINIT/m_Ctqmc/Ctqmc_getE [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_getE

FUNCTION

  get interaction energy and noise on it

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

  E=interaction energy
  noise=noise on this value

SIDE EFFECTS

NOTES

SOURCE

2760 SUBROUTINE Ctqmc_getE(this,E,noise)
2761 
2762 !Arguments ------------------------------------
2763   TYPE(Ctqmc)       , INTENT(IN ) :: this
2764   DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: E
2765   DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: Noise
2766 
2767   IF ( PRESENT(E) ) &
2768     E = this%measDE(1,1)  
2769   IF ( PRESENT(Noise) ) &
2770     Noise = SUM(this%Impurity%mat_U)/(this%flavors*(this%flavors-1)) &
2771             * this%a_Noise*(DBLE(this%sweeps)*DBLE(this%size))**this%b_Noise
2772 END SUBROUTINE Ctqmc_getE

ABINIT/m_Ctqmc/Ctqmc_getGreen [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_getGreen

FUNCTION

  Get the full green functions in time and/or frequency

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

  Gtau=green function in time
  Gw=green function in frequency

SIDE EFFECTS

NOTES

SOURCE

2629 SUBROUTINE Ctqmc_getGreen(this, Gtau, Gw)
2630 
2631 !Arguments ------------------------------------
2632   TYPE(Ctqmc)          , INTENT(INOUT)    :: this
2633   DOUBLE PRECISION, DIMENSION(:,:), OPTIONAL, INTENT(INOUT) :: Gtau
2634   COMPLEX(KIND=8), DIMENSION(:,:), OPTIONAL, INTENT(INOUT) :: Gw
2635 !Local variables ------------------------------
2636   !INTEGER                            :: itime
2637   INTEGER                            :: iflavor1
2638   INTEGER                            :: iflavor2
2639   INTEGER                            :: iflavor3
2640   INTEGER                            :: flavors
2641   DOUBLE PRECISION :: u1 
2642   DOUBLE PRECISION :: u2
2643   DOUBLE PRECISION :: Un
2644   DOUBLE PRECISION :: UUnn
2645 
2646   flavors = this%flavors
2647   DO iflavor1 = 1, flavors
2648     u1 = 0.d0
2649     u2 = 0.d0
2650     DO iflavor2 = 1, flavors
2651       IF ( iflavor2 .EQ. iflavor1 ) CYCLE
2652       Un = this%Impurity%mat_U(iflavor2,iflavor1) * this%measN(1,iflavor2)
2653       u1 = u1 + Un 
2654       u2 = u2 + Un*this%Impurity%mat_U(iflavor2,iflavor1) 
2655       DO iflavor3 = 1, flavors
2656         IF ( iflavor3 .EQ. iflavor2 .OR. iflavor3 .EQ. iflavor1 ) CYCLE
2657         UUnn = (this%Impurity%mat_U(iflavor2,iflavor1)*this%Impurity%mat_U(iflavor3,iflavor1)) * this%measDE(iflavor2,iflavor3) 
2658         u2 = u2 + UUnn 
2659       END DO
2660     END DO  
2661 
2662     CALL GreenHyb_setMoments(this%Greens(iflavor1),u1,u2)
2663     IF ( PRESENT( Gtau ) ) THEN
2664       Gtau(1:this%samples,iflavor1) = this%Greens(iflavor1)%oper(1:this%samples)
2665     END IF
2666        !write(6,*) "present gw", present(gw)
2667     IF ( PRESENT( Gw ) ) THEN
2668        !write(6,*) "size gw",SIZE(Gw,DIM=2) ,flavors+1 
2669       IF ( SIZE(Gw,DIM=2) .EQ. flavors+1 ) THEN
2670         CALL GreenHyb_forFourier(this%Greens(iflavor1), Gomega=Gw(:,iflavor1), omega=Gw(:,this%flavors+1))
2671         !IF ( this%rank .EQ. 0 ) write(20,*) Gw(:,iflavor1)
2672       ELSE IF ( SIZE(Gw,DIM=2) .EQ. flavors ) THEN  
2673         CALL GreenHyb_forFourier(this%Greens(iflavor1),Gomega=Gw(:,iflavor1))
2674       ELSE
2675         CALL WARNALL("Ctqmc_getGreen : Gw is not valid                    ")
2676         CALL GreenHyb_forFourier(this%Greens(iflavor1),Wmax=this%Wmax)
2677       END IF
2678     ELSE
2679       CALL GreenHyb_forFourier(this%Greens(iflavor1),Wmax=this%Wmax)
2680     END IF
2681   END DO
2682 END SUBROUTINE Ctqmc_getGreen

ABINIT/m_Ctqmc/Ctqmc_getResult [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_getResult

FUNCTION

  reduce everything to get the result of the simulation

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

2105 SUBROUTINE Ctqmc_getResult(this)
2106 
2107 
2108 #ifdef HAVE_MPI1
2109 include 'mpif.h'
2110 #endif
2111 !Arguments ------------------------------------
2112   TYPE(Ctqmc)  , INTENT(INOUT)                    :: this
2113 !Local variables ------------------------------
2114   INTEGER                                       :: iflavor
2115   INTEGER                                       :: flavors
2116   INTEGER                                       :: itau
2117   INTEGER                                       :: endDensity
2118   DOUBLE PRECISION                              :: inv_flavors
2119   DOUBLE PRECISION                              :: a
2120   DOUBLE PRECISION                              :: b
2121   DOUBLE PRECISION                              :: r
2122   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: alpha
2123   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: beta
2124   DOUBLE PRECISION,              DIMENSION(1:2) :: TabX
2125   DOUBLE PRECISION,              DIMENSION(1:2) :: TabY
2126   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:)   :: freqs
2127   INTEGER, ALLOCATABLE, DIMENSION(:)   :: counts
2128   INTEGER, ALLOCATABLE, DIMENSION(:)   :: displs
2129   INTEGER                                       :: sp1
2130   INTEGER                                       :: spAll
2131   INTEGER                                       :: last
2132   INTEGER                                       :: n1
2133   INTEGER                                       :: n2
2134   INTEGER                                       :: debut
2135 !  INTEGER                                       :: fin
2136 #ifdef HAVE_MPI
2137   INTEGER                                       :: ierr
2138 #endif
2139   DOUBLE PRECISION                              :: inv_size,sumh
2140   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: buffer 
2141   TYPE(FFTHyb) :: FFTmrka
2142 
2143   IF ( .NOT. this%done ) &
2144     CALL ERROR("Ctqmc_getResult : Simulation not run                ")
2145 
2146   flavors     =  this%flavors
2147   inv_flavors = 1.d0 / DBLE(flavors)
2148 
2149 
2150   inv_size = 1.d0 / DBLE(this%size)
2151   sp1 = 0
2152   spAll = 0
2153 
2154 !#ifdef CTCtqmc_CHECK
2155   IF ( this%opt_check .GT. 0 ) THEN
2156     this%errorImpurity = ImpurityOperator_getError(this%Impurity) * inv_flavors 
2157     this%errorBath     = BathOperator_getError    (this%Bath    ) * inv_flavors 
2158   END IF
2159 !#endif
2160 
2161   MALLOC(alpha,(1,1))
2162   MALLOC(beta,(1,1))
2163   MALLOC(buffer,(1,1))
2164   IF ( this%opt_noise .EQ. 1) THEN
2165     FREEIF(alpha)
2166     MALLOC(alpha,(1:this%samples+1,1:flavors))
2167     FREEIF(beta)
2168     MALLOC(beta,(1:this%samples+1,1:flavors))
2169   END IF
2170 
2171   IF ( this%have_MPI .EQV. .TRUE.) THEN 
2172     sp1   = this%samples+1
2173     spALL = sp1 + flavors + 6 
2174 
2175 !#ifdef CTCtqmc_ANALYSIS
2176     IF ( this%opt_analysis .EQ. 1 ) &
2177       spAll = spAll + 3*sp1 
2178     IF ( this%opt_order .GT. 0 ) &
2179       spAll = spAll + this%opt_order 
2180     IF ( this%opt_noise .EQ. 1 ) &
2181       spAll = spAll + 2*(this%samples + 1)
2182 !#endif
2183 
2184     FREEIF(buffer)
2185     MALLOC(buffer,(1:spAll,1:MAX(2,flavors)))
2186   END IF
2187 
2188 !  this%seg_added    = this%seg_added    * inv_flavors 
2189 !  this%seg_removed  = this%seg_removed  * inv_flavors
2190 !  this%seg_sign     = this%seg_sign     * inv_flavors
2191 !  this%anti_added   = this%anti_added   * inv_flavors
2192 !  this%anti_removed = this%anti_removed * inv_flavors
2193 !  this%anti_sign    = this%anti_sign    * inv_flavors
2194   this%stats(:) = this%stats(:) * inv_flavors
2195 
2196   DO iflavor = 1, flavors
2197     CALL GreenHyb_measHybrid(this%Greens(iflavor), this%Bath%M(iflavor), this%Impurity%Particles(iflavor), .TRUE.)
2198     CALL GreenHyb_getHybrid(this%Greens(iflavor))
2199     ! Accumule les dernieres mesure de N
2200     this%measN(1,iflavor) = this%measN(1,iflavor) + this%measN(3,iflavor)*this%measN(4,iflavor)
2201     this%measN(2,iflavor) = this%measN(2,iflavor) + this%measN(4,iflavor)
2202     ! Reduction
2203     this%measN(1,iflavor)  = this%measN(1,iflavor) / ( this%measN(2,iflavor) * this%beta )
2204     ! Correction
2205     CALL GreenHyb_setN(this%Greens(iflavor), this%measN(1,iflavor))
2206 !#ifdef CTCtqmc_ANALYSIS
2207     IF ( this%opt_order .GT. 0 ) &
2208       this%measPerturbation(:   ,iflavor) = this%measPerturbation(:,iflavor) &
2209                                     / SUM(this%measPerturbation(:,iflavor))
2210     IF ( this%opt_analysis .EQ. 1 ) THEN
2211       this%measCorrelation (:,1,iflavor) = this%measCorrelation  (:,1,iflavor) &
2212                                     / SUM(this%measCorrelation (:,1,iflavor)) &
2213                                     * this%inv_dt 
2214       this%measCorrelation (:,2,iflavor) = this%measCorrelation  (:,2,iflavor) &
2215                                     / SUM(this%measCorrelation (:,2,iflavor)) &
2216                                     * this%inv_dt 
2217       this%measCorrelation (:,3,iflavor) = this%measCorrelation  (:,3,iflavor) &
2218                                     / SUM(this%measCorrelation (:,3,iflavor)) &
2219                                     * this%inv_dt 
2220     END IF
2221 !#endif
2222     IF ( this%opt_noise .EQ. 1 ) THEN
2223       TabX(1) = DBLE(this%modNoise2)
2224       TabX(2) = DBLE(this%modNoise1)
2225       DO itau = 1, this%samples+1
2226         this%measNoiseG(itau,iflavor,2)%vec = -this%measNoiseG(itau,iflavor,2)%vec*this%inv_dt &  
2227                                            /(this%beta*DBLE(this%modNoise2))
2228         this%measNoiseG(itau,iflavor,1)%vec = -this%measNoiseG(itau,iflavor,1)%vec*this%inv_dt &  
2229                                            /(this%beta*DBLE(this%modNoise1))
2230         n2 = this%measNoiseG(itau,iflavor,2)%tail
2231         TabY(1) = Stat_deviation(this%measNoiseG(itau,iflavor,2)%vec(1:n2))!*SQRT(n2/(n2-1))
2232         n1 = this%measNoiseG(itau,iflavor,1)%tail
2233         TabY(2) = Stat_deviation(this%measNoiseG(itau,iflavor,1)%vec(1:n1))!*SQRT(n1/(n1-1))
2234         CALL Stat_powerReg(TabX,SQRT(2.d0*LOG(2.d0))*TabY,alpha(itau,iflavor),beta(itau,iflavor),r)
2235         ! ecart type -> 60%
2236         ! largeur a mi-hauteur d'une gaussienne -> sqrt(2*ln(2))*sigma
2237       END DO
2238     END IF
2239 
2240     IF ( this%have_MPI .EQV. .TRUE. ) THEN 
2241       buffer(1:sp1, iflavor) = this%Greens(iflavor)%oper(1:sp1)
2242     END IF
2243   END DO
2244   last = sp1
2245 
2246 
2247   this%measDE(:,:) = this%measDE(:,:) * DBLE(this%measurements) /(DBLE(this%sweeps)*this%beta)
2248   IF ( this%opt_histo .GT. 0 ) THEN
2249     this%occup_histo_time(:) = this%occup_histo_time(:) / INT(this%sweeps/this%measurements)
2250   END IF
2251  ! write(6,*) "=== Histogram of occupations for complete simulation ====",INT(this%sweeps/this%measurements)
2252  ! sumh=0
2253  ! do n1=1,this%flavors+1
2254  !    write(6,'(i4,f10.4)')  n1-1, this%occup_histo_time(n1)
2255  !    sumh=sumh+this%occup_histo_time(n1)
2256  ! enddo
2257  ! write(6,*) "=================================",sumh
2258 
2259   n1 = this%measNoise(1)%tail
2260   n2 = this%measNoise(2)%tail
2261 
2262   ! On utilise freqs comme tableau de regroupement
2263   ! Gather de Noise1
2264   IF ( this%have_MPI .EQV. .TRUE. ) THEN
2265     MALLOC(counts,(1:this%size))
2266     MALLOC(displs,(1:this%size))
2267     FREEIF(freqs)
2268     MALLOC(freqs,(1:this%size*n1))
2269     freqs = 0.d0
2270     freqs(n1*this%rank+1:n1*(this%rank+1)) = this%measNoise(1)%vec(1:n1) 
2271     counts(:) = n1
2272     displs(:) = (/ ( iflavor*n1, iflavor=0, this%size-1 ) /)
2273 #ifdef HAVE_MPI
2274     CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DOUBLE_PRECISION, &
2275                         freqs, counts, displs, &
2276                         MPI_DOUBLE_PRECISION, this%MY_COMM, ierr)
2277 #endif
2278     n1 = this%size*n1
2279     CALL Vector_setSize(this%measNoise(1),n1)
2280     this%measNoise(1)%vec(1:n1) = freqs(:)
2281     ! Gather de Noise2
2282     FREE(freqs)
2283     MALLOC(freqs,(1:this%size*n2))
2284     freqs = 0.d0
2285     freqs(n2*this%rank+1:n2*(this%rank+1)) = this%measNoise(2)%vec(1:n2) 
2286     counts(:) = n2
2287     displs(:) = (/ ( iflavor*n2, iflavor=0, this%size-1 ) /)
2288 #ifdef HAVE_MPI
2289     CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DATATYPE_NULL, &
2290                         freqs, counts, displs, &
2291                         MPI_DOUBLE_PRECISION, this%MY_COMM, ierr)
2292 #endif
2293     n2 = this%size*n2
2294     CALL Vector_setSize(this%measNoise(2),n2)
2295     this%measNoise(2)%vec(1:n2) = freqs(:)
2296     FREE(counts)
2297     FREE(displs)
2298     FREE(freqs)
2299   END IF
2300   !n1 = this%measNoise(1)%tail
2301   !n2 = this%measNoise(2)%tail
2302 
2303   ! Transformation des paquets pour que ca fit a CTQMC_SLICE(1|2)
2304   IF ( n1 .GT. CTQMC_SLICE1 ) THEN
2305     itau = n1/CTQMC_SLICE1
2306     MALLOC(freqs,(1:n1/itau))
2307     DO debut=1, n1/itau
2308       freqs(debut)=SUM(this%measNoise(1)%vec((debut-1)*itau+1:itau*debut))
2309     END DO
2310     freqs(:) = freqs(:)/DBLE(itau)
2311     this%modNoise1 = this%modNoise1*itau
2312     n1 = n1/itau
2313     CALL Vector_setSize(this%measNoise(1),n1)
2314     this%measNoise(1)%vec(1:n1) = freqs(:)
2315     FREE(freqs)
2316   END IF
2317   IF ( n2 .GT. CTQMC_SLICE1*CTQMC_SLICE2 ) THEN
2318     itau = n2/(CTQMC_SLICE1*CTQMC_SLICE2)
2319     MALLOC(freqs,(1:n2/itau))
2320     DO debut=1, n2/itau
2321       freqs(debut)=SUM(this%measNoise(2)%vec((debut-1)*itau+1:itau*debut))
2322     END DO
2323     freqs(:) = freqs(:)/DBLE(itau)
2324     this%modNoise2 = this%modNoise2*itau
2325     n2 = n2/itau
2326     CALL Vector_setSize(this%measNoise(2),n2)
2327     this%measNoise(2)%vec(1:n2) = freqs(:)
2328     FREE(freqs)
2329   END IF
2330   ! On peut s'amuser avec nos valeur d'energies
2331   !MALLOC(TabX,(1:20))
2332   !MALLOC(TabY,(1:20))
2333 
2334   TabX(1) = DBLE(this%modNoise2)
2335   TabX(2) = DBLE(this%modNoise1)
2336 
2337   ! Il faut calculer pour chaque modulo 10 ecarts type sur les donnes acquises
2338   this%measNoise(1)%vec(1:n1) = this%measNoise(1)%vec(1:n1)/(this%beta*DBLE(this%modNoise1))*DBLE(this%measurements)
2339   this%measNoise(2)%vec(1:n2) = this%measNoise(2)%vec(1:n2)/(this%beta*DBLE(this%modNoise2))*DBLE(this%measurements)
2340 !  CALL Vector_print(this%measNoise(1),this%rank+70)
2341 !  CALL Vector_print(this%measNoise(2),this%rank+50)
2342 !  DO iflavor=1,10
2343 !    debut = (iflavor-1)*n2/10+1
2344 !    fin   = iflavor*n2/10
2345 !    TabY(iflavor) = Stat_deviation(this%measNoise(2)%vec(debut:fin))
2346 !    debut = (iflavor-1)*n1/10+1
2347 !    fin   = iflavor*n1/10
2348 !    TabY(10+iflavor) = Stat_deviation(this%measNoise(1)%vec(debut:fin))
2349 !  END DO
2350 !!  TabY(1:n) = (this%measNoise(2)%vec(1:n)   &
2351 !!              )
2352 !!             !/(this%beta*DBLE(this%modNoise2))*DBLE(this%measurements) &
2353 !!             !- this%measDE(1,1))
2354 !!  TabY(this%measNoise(2)%tail+1:n+this%measNoise(2)%tail) = (this%measNoise(1)%vec(1:n)   &
2355 !!               )
2356 !!             ! /(this%beta*DBLE(this%modNoise1))*DBLE(this%measurements) &
2357 !!             ! - this%measDE(1,1))
2358 !  IF ( this%rank .EQ. 0 ) THEN
2359 !    DO iflavor=1,20
2360 !      write(45,*) TabX(iflavor), TabY(iflavor)
2361 !    END DO
2362 !  END IF
2363 !
2364 
2365 
2366   TabY(1) = Stat_deviation(this%measNoise(2)%vec(1:n2))!*SQRT(n2/(n2-1))
2367 !!  write(this%rank+10,*) TabX(2)
2368 !!  write(this%rank+40,*) TabX(1)
2369 !!  CALL Vector_print(this%measNoise(1),this%rank+10)
2370 !!  CALL Vector_print(this%measNoise(2),this%rank+40)
2371 !!  CLOSE(this%rank+10)
2372 !!  CLOSE(this%rank+40)
2373   TabY(2) = Stat_deviation(this%measNoise(1)%vec(1:n1))!*SQRT(n1/(n1-1))
2374 !!  ! Ecart carre moyen ~ ecart type mais non biaise. Serait moins precis. Aucun
2375   ! impact sur la pente, juste sur l'ordonnee a l'origine.
2376 
2377   CALL Stat_powerReg(TabX,SQRT(2.d0*LOG(2.d0))*TabY,a,b,r)
2378 !  FREE(TabX)
2379 !  FREE(TabY)
2380   ! ecart type -> 60%
2381   ! largeur a mi-hauteur d'une gaussienne -> sqrt(2*ln(2))*sigma
2382 
2383   !this%measDE(1,1) = SUM(this%measNoise(1)%vec(1:this%measNoise(1)%tail))/(DBLE(this%measNoise(1)%tail*this%modNoise1)*this%beta)
2384   !this%measDE(2:flavors,1:flavors) = this%measDE(2:flavors,1:flavors) /(DBLE(this%sweeps)*this%beta)
2385   CALL ImpurityOperator_getErrorOverlap(this%Impurity,this%measDE)
2386   ! Add the difference between true calculation and quick calculation of the
2387   ! last sweep overlap to measDE(2,2)
2388   !this%measDE = this%measDE * DBLE(this%measurements) 
2389   IF ( this%have_MPI .EQV. .TRUE. ) THEN 
2390     IF ( this%opt_analysis .EQ. 1 ) THEN
2391       buffer(last+1:last+sp1,:) = this%measCorrelation(:,1,:)
2392       last = last + sp1
2393       buffer(last+1:last+sp1,:) = this%measCorrelation(:,2,:)
2394       last = last + sp1
2395       buffer(last+1:last+sp1,:) = this%measCorrelation(:,3,:)
2396       last = last + sp1
2397     END IF
2398     IF ( this%opt_order .GT. 0 ) THEN
2399       buffer(last+1:last+this%opt_order, :) = this%measPerturbation(:,:)
2400       last = last + this%opt_order
2401     END IF
2402     IF ( this%opt_noise .EQ. 1 ) THEN
2403       buffer(last+1:last+this%samples+1,:) = alpha(:,:)
2404       last = last + this%samples + 1
2405       buffer(last+1:last+this%samples+1,:) = beta(:,:)
2406       last = last + this%samples + 1
2407     END IF
2408 !  this%measDE(2,2) = a*EXP(b*LOG(DBLE(this%sweeps*this%size)))
2409     buffer(spall-(flavors+5):spAll-6,:) = this%measDE(:,:)
2410 !    buffer(spAll  ,1) = this%seg_added   
2411 !    buffer(spAll-1,1) = this%seg_removed 
2412 !    buffer(spAll-2,1) = this%seg_sign    
2413 !    buffer(spAll  ,2) = this%anti_added  
2414 !    buffer(spAll-1,2) = this%anti_removed
2415 !    buffer(spAll-2,2) = this%anti_sign   
2416     buffer(spAll  ,1) = this%stats(1)
2417     buffer(spAll-1,1) = this%stats(2)
2418     buffer(spAll-2,1) = this%stats(3)
2419     buffer(spAll  ,2) = this%stats(4)
2420     buffer(spAll-1,2) = this%stats(5)
2421     buffer(spAll-2,2) = this%stats(6)
2422     buffer(spAll-3,1) = this%swap
2423     buffer(spAll-3,2) = DBLE(this%modGlobalMove(2))
2424     buffer(spAll-4,1) = a
2425     buffer(spAll-4,2) = b
2426 !#ifdef CTCtqmc_CHECK
2427     buffer(spAll-5,1) = this%errorImpurity
2428     buffer(spAll-5,2) = this%errorBath 
2429 !#endif
2430 
2431 #ifdef HAVE_MPI
2432     CALL MPI_ALLREDUCE(MPI_IN_PLACE, buffer, spAll*flavors, &
2433                      MPI_DOUBLE_PRECISION, MPI_SUM, this%MY_COMM, ierr)
2434     CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%runTime, 1, MPI_DOUBLE_PRECISION, MPI_MAX, &
2435              this%MY_COMM, ierr)
2436     IF ( this%opt_histo .GT. 0 ) THEN
2437       CALL MPI_ALLREDUCE(MPI_IN_PLACE, this%occup_histo_time, flavors+1, MPI_DOUBLE_PRECISION, MPI_SUM, &
2438                this%MY_COMM, ierr)
2439     END IF
2440 #endif
2441 
2442   
2443     buffer          = buffer * inv_size
2444     this%measDE(:,:)  = buffer(spall-(flavors+5):spAll-6,:)
2445 !    this%seg_added    = buffer(spAll  ,1)
2446 !    this%seg_removed  = buffer(spAll-1,1)
2447 !    this%seg_sign     = buffer(spAll-2,1)
2448 !    this%anti_added   = buffer(spAll  ,2)
2449 !    this%anti_removed = buffer(spAll-1,2)
2450 !    this%anti_sign    = buffer(spAll-2,2)
2451     this%stats(1)    = buffer(spAll  ,1)
2452     this%stats(2)    = buffer(spAll-1,1)
2453     this%stats(3)    = buffer(spAll-2,1)
2454     this%stats(4)    = buffer(spAll  ,2)
2455     this%stats(5)    = buffer(spAll-1,2)
2456     this%stats(6)    = buffer(spAll-2,2)
2457     this%swap         = buffer(spAll-3,1)
2458     this%modGlobalMove(2) = NINT(buffer(spAll-3,2))
2459     a               = buffer(spAll-4,1) 
2460     b               = buffer(spAll-4,2)
2461 !!#ifdef CTCtqmc_CHECK
2462     this%errorImpurity= buffer(spAll-5,1) 
2463     this%errorBath    = buffer(spAll-5,2)   
2464 !#endif
2465 
2466     DO iflavor = 1, flavors
2467       this%Greens(iflavor)%oper          = buffer(1:sp1          , iflavor)
2468     END DO
2469     last = sp1
2470     IF ( this%opt_analysis .EQ. 1 ) THEN
2471       this%measCorrelation(:,1,:) = buffer(last+1:last+sp1,:) 
2472       last = last + sp1
2473       this%measCorrelation(:,2,:) = buffer(last+1:last+sp1,:) 
2474       last = last + sp1
2475       this%measCorrelation(:,3,:) = buffer(last+1:last+sp1,:) 
2476       last = last + sp1
2477     END IF
2478     IF ( this%opt_order .GT. 0 ) THEN
2479       this%measPerturbation(:,:) = buffer(last+1:last+this%opt_order, :)
2480       last = last + this%opt_order
2481     END IF
2482     IF ( this%opt_noise .EQ. 1 ) THEN
2483       alpha(:,:) = buffer(last+1:last+this%samples+1,:)
2484       last = last + this%samples + 1
2485       beta(:,:) = buffer(last+1:last+this%samples+1,:)
2486       last = last + this%samples + 1
2487     END IF
2488   END IF
2489   DO iflavor = 1, flavors
2490     ! complete DE this
2491     this%measDE(iflavor, iflavor+1:flavors) = this%measDE(iflavor+1:flavors,iflavor)
2492   END DO
2493   FREE(buffer)
2494 
2495   IF ( this%opt_spectra .GE. 1 ) THEN
2496     endDensity = SIZE(this%density,2)
2497     IF ( this%density(1,endDensity) .EQ. -1.d0 ) &
2498       endDensity = endDensity - 1
2499     CALL FFTHyb_init(FFTmrka,endDensity,DBLE(this%thermalization)/DBLE(this%measurements*this%opt_spectra))
2500     ! Not very Beauty 
2501     MALLOC(freqs,(1:FFTmrka%size/2))
2502     DO iflavor = 1, flavors
2503       ! mean value is removed to supress the continue composent 
2504       CALL FFTHyb_setData(FFTmrka,this%density(iflavor,1:endDensity)/this%beta+this%Greens(iflavor)%oper(this%samples+1))
2505       CALL FFTHyb_run(FFTmrka,1)
2506       CALL FFTHyb_getData(FFTmrka,endDensity,this%density(iflavor,:),freqs)
2507     END DO
2508     this%density(flavors+1,:) = -1.d0
2509     this%density(flavors+1,1:FFTmrka%size/2) = freqs
2510     CALL FFTHyb_destroy(FFTmrka)
2511     FREE(freqs)
2512   END IF
2513 
2514   this%a_Noise = a
2515   this%b_Noise = b
2516   IF ( this%opt_noise .EQ. 1 ) THEN
2517     this%abNoiseG(1,:,:) = alpha
2518     this%abNoiseG(2,:,:) = beta
2519   END IF
2520   FREE(alpha)
2521   FREE(beta)
2522 
2523   IF ( this%opt_histo .GT. 0 ) THEN
2524     write(this%ostream,*) "=== Histogram of occupations for complete simulation  ===="
2525    ! write(6,*) "sumh over procs", sumh
2526     sumh=0
2527     do n1=1,this%flavors+1
2528        write(this%ostream,'(i4,f10.4)')  n1-1, this%occup_histo_time(n1)/float(this%size)
2529        sumh=sumh+this%occup_histo_time(n1)/float(this%size)
2530     enddo
2531        write(this%ostream,'(a,f10.4)') " all" , sumh
2532     write(this%ostream,*) "================================="
2533   ENDIF
2534 
2535 END SUBROUTINE Ctqmc_getResult

ABINIT/m_Ctqmc/Ctqmc_init [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_init

FUNCTION

  Initialize the type Ctqmc
  Allocate all the non optional variables

COPYRIGHT

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

INPUTS

  this=ctqmc
  ostream=where to write
  istream=where to read the input parameters if so
  bFile=do we read in istream ?
  MY_COMM=mpi communicator for the CTQMC
  iBuffer=input parameters if bFile is false

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

306 SUBROUTINE Ctqmc_init(this, ostream, istream, bFile, MY_COMM, iBuffer)
307 
308 
309 #ifdef HAVE_MPI1
310 include 'mpif.h'
311 #endif
312 !Arguments ------------------------------------
313   TYPE(Ctqmc), INTENT(INOUT)                      :: this
314   INTEGER  , INTENT(IN   )                      :: ostream
315   INTEGER  , INTENT(IN   )                      :: istream
316   LOGICAL  , INTENT(IN   )                      :: bFile
317   DOUBLE PRECISION, DIMENSION(1:9), OPTIONAL, INTENT(IN) :: iBuffer
318   INTEGER  , OPTIONAL, INTENT(IN   )                      :: MY_COMM
319 !Local variables ------------------------------
320 #ifdef HAVE_MPI
321   INTEGER                                       :: ierr
322 #endif
323   INTEGER                                       :: iflavor
324 #ifdef __GFORTRAN__
325 !  INTEGER                                       :: pid
326 !  CHARACTER(LEN=5)                              :: Cpid
327 !
328 #endif
329   DOUBLE PRECISION, DIMENSION(1:9)             :: buffer
330 
331   this%ostream = ostream
332   this%istream = istream
333   
334 ! --- RENICE ---
335 !#ifdef __GFORTRAN__
336 !  pid = GetPid()
337 !  WRITE(Cpid,'(I5)') pid
338 !  CALL SYSTEM('renice +19 '//TRIM(ADJUSTL(Cpid))//' > /dev/null')
339 !#endif
340 !! --- RENICE ---
341 
342   IF ( PRESENT(MY_COMM)) THEN
343 #ifdef HAVE_MPI
344     this%have_MPI = .TRUE.
345     this%MY_COMM = MY_COMM
346     CALL MPI_Comm_rank(this%MY_COMM, this%rank, ierr)
347     CALL MPI_Comm_size(this%MY_COMM, this%size, ierr)
348 #else
349     CALL WARN("Ctqmc_init : MPI is not used                                    ")
350     this%have_MPI = .FALSE.
351     this%MY_COMM = -1
352     this%rank = 0
353     this%size = 1
354 #endif
355   ELSE
356     this%have_MPI = .FALSE.
357     this%MY_COMM = -1
358     this%rank = 0
359     this%size = 1
360   END IF
361 
362   !IF ( this%rank .EQ. 0 ) THEN
363 !  WRITE(ostream,'(A20)') 'Job reniced with +19'
364     !CALL FLUSH(ostream)
365   !END IF
366 
367   IF ( bFile .EQV. .TRUE. ) THEN
368     IF ( this%rank .EQ. 0 ) THEN
369 
370       READ(istream,*) buffer(1) !iseed
371       READ(istream,*) buffer(2) !this%sweeps
372       READ(istream,*) buffer(3) !this%thermalization
373       READ(istream,*) buffer(4) !this%measurements
374       READ(istream,*) buffer(5) !this%flavors
375       READ(istream,*) buffer(6) !this%samples
376       READ(istream,*) buffer(7) !this%beta
377       READ(istream,*) buffer(8) !U
378       READ(istream,*) buffer(9) !iTech
379       !READ(istream,*) buffer(9) !Wmax
380 !#ifdef CTCtqmc_ANALYSIS
381       !READ(istream,*) buffer(10) !order
382 !#endif
383     END IF
384 
385 #ifdef HAVE_MPI
386     IF ( this%have_MPI .EQV. .TRUE. ) &
387       CALL MPI_Bcast(buffer, 9, MPI_DOUBLE_PRECISION, 0,    &
388                    this%MY_COMM, ierr)
389 #endif
390   ELSE IF ( PRESENT(iBuffer) ) THEN
391     buffer(1:9) = iBuffer(1:9)
392   ELSE
393     CALL ERROR("Ctqmc_init : No input parameters                    ")
394   END IF
395 
396   CALL Ctqmc_setParameters(this, buffer)
397 
398   CALL Ctqmc_allocateAll(this)
399 
400   DO iflavor = 1, this%flavors
401       CALL GreenHyb_init(this%Greens(iflavor),this%samples, this%beta, iTech=INT(buffer(9)),MY_COMM=this%MY_COMM)
402   END DO
403 
404 
405 !  this%seg_added    = 0.d0
406 !  this%anti_added   = 0.d0
407 !  this%seg_removed  = 0.d0
408 !  this%anti_removed = 0.d0
409 !  this%seg_sign     = 0.d0
410 !  this%anti_sign    = 0.d0
411   this%stats(:)     = 0.d0
412   this%swap         = 0.d0
413   this%runTime      = 0.d0
414 
415   CALL Vector_init(this%measNoise(1),this%sweeps/this%modNoise1)
416   CALL Vector_init(this%measNoise(2),(this%sweeps/this%modNoise1+1)*CTQMC_SLICE2)
417   !CALL Vector_init(this%measNoise(3),101)
418   !CALL Vector_init(this%measNoise(4),101)
419 
420   this%set  = this%para .AND. this%inF
421   this%done = .FALSE.
422   this%init = .TRUE.
423 
424 !#ifdef CTCtqmc_CHECK
425   this%errorImpurity = 0.d0
426   this%errorBath     = 0.d0
427 !#endif
428 END SUBROUTINE Ctqmc_init

ABINIT/m_Ctqmc/Ctqmc_loop [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_loop

FUNCTION

  Definition the main loop of the CT-QMC

COPYRIGHT

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

INPUTS

  this=ctqmc
  itotal=number of sweeps to perform : thermalization or sweeps
  ilatex=unit of file to write movie if so

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

1464 SUBROUTINE Ctqmc_loop(this,itotal,ilatex)
1465 
1466 !Arguments ------------------------------------
1467   TYPE(Ctqmc), INTENT(INOUT)         :: this
1468   INTEGER    , INTENT(IN   )         :: itotal
1469   INTEGER    , INTENT(IN   )         :: ilatex
1470 !Local variables ------------------------------
1471   LOGICAL                            :: updated 
1472   LOGICAL                            :: updated_seg
1473   LOGICAL, DIMENSION(:), ALLOCATABLE :: updated_swap
1474 
1475   INTEGER                            :: flavors
1476   INTEGER                            :: measurements
1477   INTEGER                            :: modNoise1
1478   INTEGER                            :: modNoise2
1479   INTEGER                            :: modGlobalMove
1480   INTEGER                            :: sp1
1481   INTEGER                            :: itau   
1482   INTEGER                            :: ind
1483   INTEGER                            :: endDensity
1484   INTEGER                            :: indDensity
1485   INTEGER                            :: swapUpdate1
1486   INTEGER                            :: swapUpdate2
1487   INTEGER                            :: old_percent
1488   INTEGER                            :: new_percent
1489   INTEGER                            :: ipercent
1490   INTEGER                            :: iflavor
1491   INTEGER                            :: isweep
1492 
1493   DOUBLE PRECISION                   :: cpu_time1
1494   DOUBLE PRECISION                   :: cpu_time2
1495   DOUBLE PRECISION                   :: NRJ_old1
1496   DOUBLE PRECISION                   :: NRJ_old2
1497   DOUBLE PRECISION                   :: NRJ_new
1498   DOUBLE PRECISION                   :: total
1499   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_old1
1500   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_old2
1501   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: gtmp_new
1502 
1503   CALL CPU_TIME(cpu_time1)
1504 
1505   flavors        = this%flavors
1506   measurements   = this%measurements
1507   modNoise1      = this%modNoise1
1508   modNoise2      = this%modNoise2
1509   modGlobalMove  = this%modGlobalMove(1)
1510   sp1            = this%samples+1
1511   IF ( this%opt_histo .GT. 0 ) THEN
1512     this%occup_histo_time= 0.d0
1513   END IF
1514   old_percent    = 0
1515 
1516   MALLOC(updated_swap,(1:flavors))
1517   updated_swap(:) = .FALSE.
1518 
1519   NRJ_old1  = 0.d0
1520   NRJ_old2  = 0.d0
1521   NRJ_new   = 0.d0
1522 
1523   MALLOC(gtmp_new,(1,1))
1524   gtmp_new  = 0.d0
1525   MALLOC(gtmp_old1,(1,1))
1526   gtmp_old1 = 0.d0
1527   MALLOC(gtmp_old2,(1,1))
1528   gtmp_old2 = 0.d0
1529 
1530   endDensity = SIZE(this%density,2)
1531 
1532   IF ( this%opt_noise .GT. 0 ) THEN
1533     FREEIF(gtmp_new)
1534     MALLOC(gtmp_new,(1:sp1,1:flavors))
1535     FREEIF(gtmp_old1)
1536     MALLOC(gtmp_old1,(1:sp1,1:flavors))
1537     FREEIF(gtmp_old2)
1538     MALLOC(gtmp_old2,(1:sp1,1:flavors))
1539   END IF
1540 
1541   IF ( this%rank .EQ. 0 ) THEN
1542     WRITE(this%ostream, '(1x,103A)') &
1543     "|----------------------------------------------------------------------------------------------------|"
1544     WRITE(this%ostream,'(1x,A)', ADVANCE="NO") "|"
1545   END IF
1546 
1547   total = DBLE(itotal)
1548 
1549   indDensity = 1
1550   DO isweep = 1, itotal
1551     DO iflavor = 1, flavors
1552       ImpurityOperator_QuickActivation(this%Impurity,iflavor)
1553       BathOperator_QuickActivation(this%Bath,iflavor)
1554 
1555       CALL Ctqmc_tryAddRemove(this,updated_seg)
1556       updated = updated_seg .OR.  updated_swap(iflavor)
1557       updated_swap(iflavor) = .FALSE.
1558 
1559       CALL GreenHyb_measHybrid(this%Greens(iflavor), this%Bath%M(iflavor), this%Impurity%Particles(iflavor), updated)
1560       CALL Ctqmc_measN        (this, iflavor, updated)
1561       IF ( this%opt_analysis .EQ. 1 ) &
1562         CALL Ctqmc_measCorrelation (this, iflavor)
1563       IF ( this%opt_order .GT. 0 ) &
1564         CALL Ctqmc_measPerturbation(this, iflavor)
1565     END DO
1566 
1567     IF ( MOD(isweep,modGlobalMove) .EQ. 0 ) THEN
1568       CALL Ctqmc_trySwap(this,swapUpdate1, swapUpdate2)
1569       IF ( swapUpdate1 .NE. 0 .AND. swapUpdate2 .NE. 0 ) THEN
1570         updated_swap(swapUpdate1) = .TRUE.
1571         updated_swap(swapUpdate2) = .TRUE.
1572       END IF
1573     END IF
1574     
1575     IF ( MOD(isweep,measurements) .EQ. 0 ) THEN
1576       CALL ImpurityOperator_measDE(this%Impurity,this%measDE)
1577       IF ( this%opt_spectra .GE. 1) THEN
1578         IF ( measurements*this%opt_spectra /= 0) THEN
1579           IF ( MOD(isweep,measurements*this%opt_spectra) .EQ. 0 ) THEN
1580              this%density(1:flavors,indDensity) = this%measN(3,1:flavors)
1581              indDensity = indDensity+1
1582           END IF
1583         END IF
1584       END IF
1585     END IF
1586 
1587     IF ( MOD(isweep,measurements) .EQ. 0 ) THEN
1588       IF ( this%opt_histo .GT. 0 ) THEN
1589         CALL ImpurityOperator_occup_histo_time(this%Impurity,this%occup_histo_time)
1590       ENDIF
1591     ENDIF
1592 
1593     IF ( MOD(isweep, modNoise1) .EQ. 0 ) THEN
1594       !modNext = isweep + modNoise2
1595       NRJ_new = (SUM(this%measDE(:,:))-this%measDE(1,1))*0.5d0 ! double occupation, avoid stat with 0 for U=J=0
1596       CALL Vector_pushBack(this%measNoise(1),NRJ_new - NRJ_old1)
1597       NRJ_old1 = NRJ_new
1598 
1599       !! Try to limit accumulation error
1600       CALL ImpurityOperator_cleanOverlaps(this%Impurity)
1601 
1602       IF ( this%opt_noise .EQ. 1 ) THEN
1603         DO iflavor = 1, flavors
1604           DO ind = 1, this%Greens(iflavor)%this%tail
1605             itau = this%Greens(iflavor)%this%listINT(ind)
1606             gtmp_new(itau,iflavor) = this%Greens(iflavor)%oper(itau) & 
1607                         +this%Greens(iflavor)%this%listDBLE(ind)*DBLE(this%Greens(iflavor)%factor)
1608           END DO
1609           DO itau = 1, sp1
1610             CALL Vector_pushBack(this%measNoiseG(itau,iflavor,1), gtmp_new(itau,iflavor) - gtmp_old1(itau,iflavor))
1611             gtmp_old1(itau,iflavor) = gtmp_new(itau,iflavor)
1612           END DO
1613         END DO
1614       END IF
1615     END IF
1616 
1617     IF ( MOD(isweep,modNoise2) .EQ. 0 ) THEN
1618       NRJ_new = (SUM(this%measDE(:,:))-this%measDE(1,1))*0.5d0 ! double occupation, avoid stat with 0 for U=J=0
1619       CALL Vector_pushBack(this%measNoise(2),NRJ_new - NRJ_old2)
1620       NRJ_old2 = NRJ_new
1621       IF ( this%opt_noise .EQ. 1 ) THEN
1622         DO iflavor = 1, flavors
1623           DO ind = 1, this%Greens(iflavor)%this%tail
1624             itau = this%Greens(iflavor)%this%listINT(ind)
1625             gtmp_new(itau,iflavor) = this%Greens(iflavor)%oper(itau) & 
1626                         +this%Greens(iflavor)%this%listDBLE(ind)*this%Greens(iflavor)%factor
1627           END DO
1628           DO itau = 1, sp1
1629             CALL Vector_pushBack(this%measNoiseG(itau,iflavor,2), gtmp_new(itau,iflavor) - gtmp_old2(itau,iflavor))
1630             gtmp_old2(itau,iflavor) = gtmp_new(itau,iflavor)
1631           END DO
1632         END DO 
1633       END IF
1634 
1635       IF ( this%rank .EQ. 0 ) THEN 
1636         new_percent = CEILING(DBLE(isweep)*100.d0/DBLE(itotal))
1637         DO ipercent = old_percent+1, new_percent 
1638           WRITE(this%ostream,'(A)',ADVANCE="NO") "-"
1639         END DO
1640         old_percent = new_percent
1641       END IF
1642     END IF
1643 
1644     IF ( this%opt_movie .EQ. 1 ) THEN
1645       WRITE(ilatex,'(A11,I9)') "%iteration ", isweep
1646       CALL ImpurityOperator_printLatex(this%Impurity,ilatex,isweep)
1647     END IF
1648 
1649   END DO
1650 
1651   IF ( this%rank .EQ. 0 ) THEN
1652     DO ipercent = old_percent+1, 100
1653       WRITE(this%ostream,'(A)',ADVANCE="NO") "-"
1654     END DO
1655     WRITE(this%ostream,'(A)') "|"
1656   END IF
1657  
1658   FREE(gtmp_new)
1659   FREE(gtmp_old1)
1660   FREE(gtmp_old2)
1661   FREE(updated_swap)
1662 
1663   IF ( this%opt_spectra .GE. 1 .AND. itotal .EQ. this%sweeps ) THEN
1664     IF ( endDensity .NE. indDensity-1 ) THEN
1665       this%density(:,endDensity) = -1.d0
1666     END IF
1667   END IF
1668 
1669   CALL CPU_TIME(cpu_time2)
1670 
1671   this%runTime = (cpu_time2 - cpu_time1)*1.05d0 ! facteur arbitraire de correction
1672 END SUBROUTINE Ctqmc_loop

ABINIT/m_Ctqmc/Ctqmc_measCorrelation [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_measCorrelation

FUNCTION

  measure all correlations in times for a flavor

COPYRIGHT

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

INPUTS

  this=ctqmc
  iflavor=the flavor to measure

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

1977 SUBROUTINE Ctqmc_measCorrelation(this, iflavor)
1978 
1979 !Arguments ------------------------------------
1980   TYPE(Ctqmc)             , INTENT(INOUT)       :: this
1981   !TYPE(ImpurityOperator), INTENT(IN   )       :: impurity
1982   INTEGER               , INTENT(IN   )       :: iflavor
1983 !Local variables ------------------------------
1984   INTEGER                                     :: iCdag
1985   INTEGER                                     :: iCdagBeta
1986   INTEGER                                     :: iC
1987   INTEGER                                     :: index
1988   INTEGER                                     :: size
1989   DOUBLE PRECISION                            :: tC
1990   DOUBLE PRECISION                            :: tCdag
1991   !DOUBLE PRECISION                            :: time
1992   DOUBLE PRECISION                            :: inv_dt
1993   DOUBLE PRECISION                            :: beta
1994 
1995   IF ( .NOT. this%set ) &
1996     CALL ERROR("Ctqmc_measCorrelation : QMC not set                 ")
1997 
1998   size = this%impurity%particles(this%impurity%activeFlavor)%tail
1999   beta = this%beta
2000 
2001   IF ( size .EQ. 0 ) RETURN
2002   
2003   inv_dt = this%inv_dt
2004 
2005   DO iCdag = 1, size ! first segments
2006     tCdag  = this%impurity%particles(this%impurity%activeFlavor)%list(iCdag,Cdag_)
2007     tC     = this%impurity%particles(this%impurity%activeFlavor)%list(iCdag,C_   )
2008     index = INT( ( (tC - tCdag)  * inv_dt ) + .5d0 ) + 1
2009     this%measCorrelation(index,1,iflavor) = this%measCorrelation(index,1,iflavor) + 1.d0
2010     MODCYCLE(iCdag+1,size,iCdagBeta)
2011     index = INT( ( ( &
2012                     this%impurity%particles(this%impurity%activeFlavor)%list(iCdagBeta,Cdag_) - tC &
2013                     + AINT(DBLE(iCdag)/DBLE(size))*beta &
2014                    )  * inv_dt ) + .5d0 ) + 1
2015     IF ( index .LT. 1 .OR. index .GT. this%samples+1 ) THEN
2016       CALL WARN("Ctqmc_measCorrelation : bad index line 1095         ")
2017     ELSE
2018       this%measCorrelation(index,2,iflavor) = this%measCorrelation(index,2,iflavor) + 1.d0
2019     END IF
2020 !    DO iC = 1, size
2021 !      tC = impurity%particles(impurity%activeFlavor)%list(C_,iC)
2022 !      time = tC - tCdag
2023 !      IF ( time .LT. 0.d0 ) time = time + beta
2024 !      index = INT( ( time * inv_dt ) + .5d0 ) + 1
2025 !      this%measCorrelation(index,3,iflavor) = this%measCorrelation(index,3,iflavor) + 1.d0
2026 !    END DO
2027     DO iC = 1, size!  this%Greens(iflavor)%index_old%tail 
2028         this%measCorrelation(this%Greens(iflavor)%this%listINT(iC+(iCdag-1)*size),3,iflavor) = &
2029         this%measCorrelation(this%Greens(iflavor)%this%listINT(iC+(iCdag-1)*size),3,iflavor) + 1.d0
2030     END DO
2031   END DO
2032 
2033 END SUBROUTINE Ctqmc_measCorrelation

ABINIT/m_Ctqmc/Ctqmc_measN [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_measN

FUNCTION

  measure the number of electron

COPYRIGHT

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

INPUTS

  this=ctqmc
  iflavor=which flavor to measure
  updated=something has changed since last time

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

1928 SUBROUTINE Ctqmc_measN(this, iflavor, updated)
1929 
1930 !Arguments ------------------------------------
1931   TYPE(Ctqmc)             , INTENT(INOUT)     :: this
1932   !TYPE(ImpurityOperator), INTENT(IN   )     :: impurity
1933   INTEGER               , INTENT(IN   )     :: iflavor
1934   LOGICAL               , INTENT(IN   )     :: updated
1935 
1936 !  IF ( .NOT. this%set ) &
1937 !    CALL ERROR("Ctqmc_measN : QMC not set                           ")
1938 
1939   
1940   IF ( updated .EQV. .TRUE. ) THEN
1941     this%measN(1,iflavor) = this%measN(1,iflavor) + this%measN(3,iflavor)*this%measN(4,iflavor)
1942     this%measN(2,iflavor) = this%measN(2,iflavor) + this%measN(4,iflavor)
1943     this%measN(3,iflavor) = ImpurityOperator_measN(this%impurity)
1944     this%measN(4,iflavor) = 1.d0
1945   ELSE
1946     this%measN(4,iflavor) = this%measN(4,iflavor) + 1.d0
1947   END IF
1948 END SUBROUTINE Ctqmc_measN

ABINIT/m_Ctqmc/Ctqmc_measPerturbation [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_measPerturbation

FUNCTION

  measure perturbation order

COPYRIGHT

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

INPUTS

  this=ctqmc
  iflavor=the flavor to measure

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

2061 SUBROUTINE Ctqmc_measPerturbation(this, iflavor)
2062 
2063 !Arguments ------------------------------------
2064   TYPE(Ctqmc)             , INTENT(INOUT)     :: this
2065   !TYPE(ImpurityOperator), INTENT(IN   )     :: impurity
2066   INTEGER               , INTENT(IN   )     :: iflavor
2067 !Local variables ------------------------------
2068   INTEGER                                   :: index
2069 
2070   IF ( .NOT. this%set ) &
2071     CALL ERROR("Ctqmc_measiPerturbation : QMC not set               ")
2072 
2073   index = this%impurity%particles(this%impurity%activeFlavor)%tail + 1
2074   IF ( index .LE. this%opt_order ) &
2075     this%measPerturbation(index,iflavor) = this%measPerturbation(index,iflavor) + 1.d0
2076 
2077 END SUBROUTINE Ctqmc_measPerturbation

ABINIT/m_Ctqmc/Ctqmc_printAll [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_printAll

FUNCTION

  print different functions computed during the simulation

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

2799 SUBROUTINE Ctqmc_printAll(this)
2800 
2801 !Arguments ------------------------------------
2802   TYPE(Ctqmc), INTENT(INOUT) :: this
2803 
2804   IF ( .NOT. this%done ) &
2805     CALL WARNALL("Ctqmc_printAll : Simulation not run                 ")
2806 
2807   CALL Ctqmc_printQMC(this)
2808 
2809   CALL Ctqmc_printGreen(this)
2810 
2811   CALL Ctqmc_printD(this)
2812 
2813 !  CALL Ctqmc_printE(this)
2814 
2815 !#ifdef CTCtqmc_ANALYSIS
2816   CALL Ctqmc_printPerturbation(this)
2817 
2818   CALL Ctqmc_printCorrelation(this)
2819 !#endif
2820 
2821   CALL Ctqmc_printSpectra(this)
2822 
2823 END SUBROUTINE Ctqmc_printAll

ABINIT/m_Ctqmc/Ctqmc_printCorrelation [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_printCorrelation

FUNCTION

  print correlation fonctions

COPYRIGHT

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

INPUTS

  this=ctqmc
  oFileIn=file stream

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

3207 SUBROUTINE Ctqmc_printCorrelation(this, oFileIn)
3208 
3209 !Arguments ------------------------------------
3210   TYPE(Ctqmc)          , INTENT(IN)             :: this
3211   INTEGER  , OPTIONAL, INTENT(IN)             :: oFileIn
3212 !Local variables ------------------------------
3213   INTEGER                                     :: oFile
3214   INTEGER                                     :: itime
3215   INTEGER                                     :: sp1
3216   INTEGER                                     :: iflavor
3217   INTEGER                                     :: i
3218   INTEGER                                     :: flavors
3219   CHARACTER(LEN=2)                            :: a
3220   CHARACTER(LEN=50)                           :: string
3221   DOUBLE PRECISION                            :: dt
3222 
3223   !IF ( this%rank .NE. MOD(5,this%size)) RETURN
3224   IF ( this%rank .NE. MOD(this%size+5,this%size)) RETURN
3225   IF ( this%opt_analysis .NE. 1 ) RETURN
3226 
3227   oFile = 44
3228   IF ( PRESENT(oFileIn) ) THEN
3229     oFile = oFileIn
3230   ELSE
3231     OPEN(UNIT=oFile, FILE="Correlation.dat")
3232   END IF
3233 
3234   sp1         =  this%samples
3235   dt          =  this%beta / sp1
3236   sp1         =  sp1 + 1
3237   flavors     =  this%flavors
3238 
3239   i = 3*flavors + 1
3240   WRITE(a,'(I2)') i
3241   WRITE(oFile,*) "# time  (/ (segement, antiseg, correl), i=1, flavor/)"
3242   string = '(1x,'//TRIM(ADJUSTL(a))//'F19.15)'
3243   DO itime = 1, sp1
3244     WRITE(oFile,string) DBLE(itime-1)*dt, &
3245                    (/ ( &
3246                    (/ ( this%measCorrelation(itime, i, iflavor), i=1,3) /) &
3247                    , iflavor=1, flavors) /)
3248   END DO
3249 
3250   IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
3251 
3252 END SUBROUTINE Ctqmc_printCorrelation

ABINIT/m_Ctqmc/Ctqmc_printD [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_printD

FUNCTION

  print individual double occupancy

COPYRIGHT

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

INPUTS

  this=ctqmc
  oFileIn=file stream

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

3031 SUBROUTINE Ctqmc_printD(this,oFileIn)
3032 
3033 !Arguments ------------------------------------
3034   TYPE(Ctqmc)          , INTENT(IN)    :: this
3035   INTEGER  , OPTIONAL, INTENT(IN)    :: oFileIn
3036 !Local variables ------------------------------
3037   INTEGER                            :: oFile
3038   INTEGER                            :: iflavor1
3039   INTEGER                            :: iflavor2
3040 
3041   !IF ( this%rank .NE. MOD(2,this%size)) RETURN
3042   IF ( this%rank .NE. MOD(this%size+2,this%size)) RETURN
3043 
3044   oFile = 41
3045   IF ( PRESENT(oFileIn) ) THEN
3046     oFile = oFileIn
3047   ELSE
3048     OPEN(UNIT=oFile, FILE="D.dat")
3049   END IF
3050 
3051   DO iflavor1 = 1, this%flavors
3052     DO iflavor2 = iflavor1+1, this%flavors
3053       WRITE(oFile,'(1x,A8,I4,A1,I4,A3,ES21.14)') "Orbitals", iflavor1, "-", iflavor2, " : ", this%measDE(iflavor2,iflavor1)
3054     END DO
3055   END DO
3056 
3057   IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
3058 
3059 END SUBROUTINE Ctqmc_printD

ABINIT/m_Ctqmc/Ctqmc_printE [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_printE

FUNCTION

  print energy and noise 

COPYRIGHT

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

INPUTS

  this=ctqmc
  oFileIn=file stream

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

3087 SUBROUTINE Ctqmc_printE(this,oFileIn)
3088 
3089 !Arguments ------------------------------------
3090   TYPE(Ctqmc)          , INTENT(IN)    :: this
3091   INTEGER  , OPTIONAL, INTENT(IN)    :: oFileIn
3092 !Local variables ------------------------------
3093   INTEGER                            :: oFile
3094   DOUBLE PRECISION                   :: E
3095   DOUBLE PRECISION                   :: Noise
3096 
3097   !IF ( this%rank .NE. MOD(3,this%size)) RETURN
3098   IF ( this%rank .NE. MOD(this%size+3,this%size)) RETURN
3099 
3100   oFile = 42
3101   IF ( PRESENT(oFileIn) ) THEN
3102     oFile = oFileIn
3103   ELSE
3104     OPEN(UNIT=oFile, FILE="BetaENoise.dat")
3105   END IF
3106 
3107   CALL Ctqmc_getE(this,E,Noise)
3108 
3109   WRITE(oFile,'(1x,F5.2,A2,ES21.14,A2,ES21.14)') this%beta, "  ", E, "  ",  Noise
3110 
3111   IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
3112 
3113 END SUBROUTINE Ctqmc_printE

ABINIT/m_Ctqmc/Ctqmc_printGreen [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_printGreen

FUNCTION

  print green functions

COPYRIGHT

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

INPUTS

  this=ctqmc
  oFileIn=file stream

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

2952 SUBROUTINE Ctqmc_printGreen(this, oFileIn)
2953 
2954 !Arguments ------------------------------------
2955   TYPE(Ctqmc)        , INTENT(IN)    :: this
2956   INTEGER  , OPTIONAL, INTENT(IN)    :: oFileIn
2957 !Local variables ------------------------------
2958   INTEGER                            :: oFile
2959   INTEGER                            :: itime
2960   INTEGER                            :: sp1
2961   INTEGER                            :: iflavor
2962   INTEGER                            :: flavors
2963   CHARACTER(LEN=4)                   :: cflavors
2964   CHARACTER(LEN=50)                  :: string
2965   DOUBLE PRECISION                   :: dt
2966   DOUBLE PRECISION                   :: sweeps
2967 
2968   !IF ( this%rank .NE. MOD(1,this%size)) RETURN
2969   IF ( this%rank .NE. MOD(this%size+1,this%size)) RETURN
2970 
2971   oFile = 40
2972   IF ( PRESENT(oFileIn) ) THEN
2973     oFile = oFileIn
2974   ELSE
2975     OPEN(UNIT=oFile, FILE="Gtau.dat")
2976   END IF
2977 
2978   sp1     =  this%samples
2979   dt      =  this%beta / DBLE(sp1)
2980   sp1     =  sp1 + 1
2981   flavors =  this%flavors
2982   sweeps = DBLE(this%sweeps)*DBLE(this%size)
2983 
2984   IF ( this%opt_noise .EQ. 1) THEN
2985     WRITE(cflavors,'(I4)') 2*flavors+1
2986     string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)'
2987     DO itime = 1, sp1
2988       WRITE(oFile,string) DBLE(itime-1)*dt, &
2989       (/ (this%Greens(iflavor)%oper(itime), iflavor=1, flavors) /), &
2990       (/ (this%abNoiseG(1,itime,iflavor)*(sweeps)**this%abNoiseG(2,itime,iflavor), iflavor=1, flavors) /)
2991     END DO
2992   ELSE
2993     WRITE(cflavors,'(I4)') flavors+1
2994     string = '(1x,'//TRIM(ADJUSTL(cflavors))//'ES22.14)'
2995     DO itime = 1, sp1
2996       WRITE(oFile,string) DBLE(itime-1)*dt, &
2997       (/ (this%Greens(iflavor)%oper(itime), iflavor=1, flavors) /)
2998     END DO
2999   END IF
3000 
3001   IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
3002 
3003 END SUBROUTINE Ctqmc_printGreen

ABINIT/m_Ctqmc/Ctqmc_printPerturbation [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_printPerturbation

FUNCTION

  print perturbation order

COPYRIGHT

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

INPUTS

  this=ctqmc
  oFileIn=file stream

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

SOURCE

3143 SUBROUTINE Ctqmc_printPerturbation(this, oFileIn)
3144 
3145 !Arguments ------------------------------------
3146   TYPE(Ctqmc)          , INTENT(IN)           :: this
3147   INTEGER  , OPTIONAL,  INTENT(IN)          :: oFileIn
3148 !Local variables-------------------------------
3149   INTEGER                                   :: oFile
3150   INTEGER                                   :: iorder
3151   INTEGER                                   :: order
3152   INTEGER                                   :: iflavor
3153   INTEGER                                   :: flavors
3154   CHARACTER(LEN=2)                          :: a
3155   CHARACTER(LEN=50)                         :: string
3156 
3157   !IF ( this%rank .NE. MOD(4,this%size)) RETURN
3158   IF ( this%rank .NE. MOD(this%size+4,this%size)) RETURN
3159   IF ( this%opt_order .LE. 0 ) RETURN
3160 
3161   oFile = 43
3162   IF ( PRESENT(oFileIn) ) THEN
3163     oFile = oFileIn
3164   ELSE
3165     OPEN(UNIT=oFile, FILE="Perturbation.dat")
3166   END IF
3167     
3168   order        =  this%opt_order
3169   flavors      =  this%flavors
3170 
3171   WRITE(a,'(I2)') flavors
3172   string = '(I5,'//TRIM(ADJUSTL(a))//'F19.15)'
3173   DO iorder = 1, order
3174     WRITE(oFile,string) iorder-1, &
3175                 (/ (this%measPerturbation(iorder, iflavor), iflavor=1, flavors) /)
3176   END DO
3177 
3178   IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
3179 END SUBROUTINE Ctqmc_printPerturbation

ABINIT/m_Ctqmc/Ctqmc_printQMC [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_printQMC

FUNCTION

  print ctqmc statistics

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

2850 SUBROUTINE Ctqmc_printQMC(this)
2851 
2852 !Arguments ------------------------------------
2853   TYPE(Ctqmc), INTENT(INOUT) :: this
2854 !Local variables ------------------------------
2855   INTEGER                  :: ostream
2856   INTEGER                  :: iflavor
2857   DOUBLE PRECISION         :: sweeps
2858   DOUBLE PRECISION         :: invSweeps
2859   CHARACTER(LEN=2)         :: a
2860   CHARACTER(LEN=15)        :: string
2861 
2862   !IF ( this%rank .NE. 0) RETURN
2863   IF ( this%rank .NE. MOD(this%size,this%size)) RETURN
2864 
2865   ostream   = this%ostream
2866   sweeps    = DBLE(this%sweeps)
2867   invSweeps = 1.d0/sweeps
2868 
2869   WRITE(ostream,'(1x,F13.0,A11,F10.2,A12,I5,A5)') sweeps*DBLE(this%size), " sweeps in ", this%runTime, &
2870                  " seconds on ", this%size, " CPUs"
2871   WRITE(ostream,'(A28,F6.2)') "Segments added        [%] : ", this%stats(CTQMC_SEGME+CTQMC_ADDED)*invSweeps*100.d0
2872   WRITE(ostream,'(A28,F6.2)') "Segments removed      [%] : ", this%stats(CTQMC_SEGME+CTQMC_REMOV)*invSweeps*100.d0
2873   WRITE(ostream,'(A28,F6.2)') "Segments sign         [%] : ", this%stats(CTQMC_SEGME+CTQMC_DETSI)*invSweeps*100.d0
2874   WRITE(ostream,'(A28,F6.2)') "Anti-segments added   [%] : ", this%stats(CTQMC_ANTIS+CTQMC_ADDED)*invSweeps*100.d0
2875   WRITE(ostream,'(A28,F6.2)') "Anti-segments removed [%] : ", this%stats(CTQMC_ANTIS+CTQMC_REMOV)*invSweeps*100.d0
2876   WRITE(ostream,'(A28,F6.2)') "Anti-segments sign    [%] : ", this%stats(CTQMC_ANTIS+CTQMC_DETSI)*invSweeps*100.d0
2877   IF ( this%modGlobalMove(1) .LT. this%sweeps + 1 ) THEN
2878     WRITE(ostream,'(A28,F6.2)') "Global Move           [%] : ", this%swap         *invSweeps*100.d0*this%modGlobalMove(1)
2879     WRITE(ostream,'(A28,F6.2)') "Global Move Reduced   [%] : ", this%swap         / DBLE(this%modGlobalMove(2))*100.d0
2880   END IF
2881 !#ifdef CTCtqmc_CHECK
2882   IF ( this%opt_check .EQ. 1 .OR. this%opt_check .EQ. 3 ) &
2883     WRITE(ostream,'(A28,E22.14)') "Impurity test         [%] : ", this%errorImpurity*100.d0
2884   IF ( this%opt_check .GE. 2 ) &
2885       WRITE(ostream,'(A28,E22.14)') "Bath     test         [%] : ", this%errorBath    *100.d0
2886 !#endif
2887   WRITE(ostream,'(A28,ES22.14,A5,ES21.14)') "<Epot>                [U] : ", this%measDE(1,1), " +/- ",&
2888 !#ifdef HAVE_MPI
2889          SUM(this%Impurity%mat_U)/(this%flavors*(this%flavors-1)) * this%a_Noise*(sweeps*DBLE(this%size))**this%b_Noise
2890 !#else
2891 !                                                              this%a_Noise*(sweeps)**this%b_Noise
2892 !#endif
2893   WRITE(ostream,'(A28,F8.4,A3,F7.4)') "Noise                [/U] : ", this%a_Noise, " x^", this%b_Noise
2894   WRITE(ostream,'(A28,E10.2)')  "Niquist puls.     [/beta] : ", ACOS(-1.d0)*this%inv_dt
2895   WRITE(ostream,'(A28,E22.14)') "Max Acc. Epot Error   [U] : ", this%measDE(2,2)/(this%beta*this%modNoise1*2.d0)*sweeps
2896   
2897   !WRITE(ostream,'(A28,F7.4,A3,F7.4,A4,E20.14)') "Noise            [G(tau)] : ", this%a_Noise(2), "x^", this%b_Noise(2), " -> ", &
2898                                                               !this%a_Noise(2)*(sweeps*DBLE(this%size))**this%b_Noise(2)
2899   IF ( this%opt_order .GT. 0 ) THEN 
2900     WRITE(a,'(I2)') this%flavors
2901     string = '(A28,'//TRIM(ADJUSTL(a))//'(1x,I3))'
2902     WRITE(ostream,string) "Perturbation orders       : ", &
2903       (/ (MAXLOC(this%measPerturbation(:, iflavor))-1, iflavor=1, this%flavors) /)
2904   END IF
2905   !CALL FLUSH(this%ostream)
2906   IF ( ABS(((this%stats(CTQMC_SEGME+CTQMC_ADDED) *invSweeps*100.d0) / &
2907             (this%stats(CTQMC_SEGME+CTQMC_REMOV) *invSweeps*100.d0) - 1.d0)) .GE. 0.02d0 &
2908    .OR. ABS(((this%stats(CTQMC_ANTIS+CTQMC_ADDED)*invSweeps*100.d0) / &
2909              (this%stats(CTQMC_ANTIS+CTQMC_REMOV)*invSweeps*100.d0) - 1.d0)) .GE. 0.02d0 ) &
2910     THEN 
2911     CALL WARNALL("Ctqmc_printQMC : bad statistic according to moves. Increase sweeps")
2912   END IF
2913   ! Check sign problem for diagonal hybridization.
2914   IF ( (this%stats(CTQMC_SEGME+CTQMC_DETSI) + this%stats(CTQMC_ANTIS+CTQMC_DETSI)) .GT. 1.d-10 ) THEN
2915     CALL WARNALL("Ctqmc_printQMC : at least one negative sign occured. There might be a bug in the CT-QMC")
2916   END IF
2917 
2918   IF ( ABS(this%b_Noise+0.5)/0.5d0 .GE. 0.05d0 ) &
2919     CALL WARNALL("Ctqmc_printQMC : bad statistic according to Noise. Increase sweeps")
2920 !  IF ( ISNAN(this%a_Noise) .OR. ISNAN(this%a_Noise) ) &
2921 !    CALL WARNALL("Ctqmc_printQMC : NaN appeared. Increase sweeps    ")
2922 
2923 
2924 END SUBROUTINE Ctqmc_printQMC

ABINIT/m_Ctqmc/Ctqmc_printSpectra [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_printSpectra

FUNCTION

  print fourier transform of time evolution of number of electrons

COPYRIGHT

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

INPUTS

  this=ctqmc
  oFileIn=file stream

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

3281 SUBROUTINE Ctqmc_printSpectra(this, oFileIn)
3282 
3283 !Arguments ------------------------------------
3284   TYPE(Ctqmc)          , INTENT(IN)             :: this
3285   INTEGER  , OPTIONAL, INTENT(IN)             :: oFileIn
3286 !Local variables ------------------------------
3287   INTEGER                                     :: oFile
3288   INTEGER                                     :: flavors
3289   INTEGER                                     :: indDensity
3290   INTEGER                                     :: endDensity
3291   CHARACTER(LEN=4)                            :: a
3292   CHARACTER(LEN=16)                           :: formatSpectra
3293 
3294   !IF ( this%rank .NE. MOD(6,this%size)) RETURN
3295   IF ( this%opt_spectra .LT. 1 ) RETURN
3296 
3297   oFile = 45+this%rank
3298   a ="0000"
3299   WRITE(a,'(I4)') this%rank
3300   IF ( PRESENT(oFileIn) ) THEN
3301     oFile = oFileIn
3302   ELSE
3303     OPEN(UNIT=oFile, FILE="Markov_"//TRIM(ADJUSTL(a))//".dat")
3304   END IF
3305 
3306   flavors     =  this%flavors
3307   WRITE(a,'(I4)') flavors+1
3308   formatSpectra ='(1x,'//TRIM(ADJUSTL(a))//'ES22.14)'
3309   WRITE(oFile,*) "# freq[/hermalization] FFT"
3310 
3311   endDensity = SIZE(this%density,2)
3312   DO WHILE ( this%density(flavors+1,endDensity) .EQ. -1 )
3313     endDensity = endDensity -1
3314   END DO
3315 
3316   DO indDensity = 1, endDensity
3317     WRITE(oFile,formatSpectra) this%density(flavors+1,indDensity), this%density(1:flavors,indDensity)
3318   END DO
3319 
3320   IF ( .NOT. PRESENT(oFileIn) ) CLOSE(oFile)
3321 
3322 END SUBROUTINE Ctqmc_printSpectra

ABINIT/m_Ctqmc/Ctqmc_reset [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_reset

FUNCTION

  reset a ctqmc simulation

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

1012 SUBROUTINE Ctqmc_reset(this)
1013 
1014 !Arguments ------------------------------------
1015   TYPE(Ctqmc), INTENT(INOUT) :: this
1016 !Local variables ------------------------------
1017   INTEGER                  :: iflavor
1018   DOUBLE PRECISION         :: sweeps
1019 
1020   DO iflavor = 1, this%flavors
1021     CALL GreenHyb_reset(this%Greens(iflavor))
1022   END DO
1023   CALL Ctqmc_clear(this)
1024   CALL ImpurityOperator_reset(this%Impurity)
1025   CALL BathOperator_reset    (this%Bath)
1026   this%measN(3,:) = 0.d0
1027   !complete restart -> measN=0
1028   this%done = .FALSE.
1029   this%set  = .FALSE.
1030   this%inF  = .FALSE.
1031   this%opt_movie = 0
1032   this%opt_analysis = 0
1033   this%opt_order = 0
1034   this%opt_check = 0
1035   this%opt_noise = 0
1036   this%opt_spectra = 0
1037   this%opt_levels = 0
1038   sweeps = DBLE(this%sweeps)*DBLE(this%size)
1039   CALL Ctqmc_setSweeps(this, sweeps)
1040 !#ifdef HAVE_MPI
1041 !  CALL MPI_BARRIER(this%MY_COMM,iflavor)
1042 !  IF ( this%rank .EQ. 0 ) &
1043 !#endif
1044 !  WRITE(this%ostream,'(A9)') "QMC reset"
1045 !  CALL FLUSH(this%ostream)
1046 END SUBROUTINE Ctqmc_reset

ABINIT/m_Ctqmc/Ctqmc_run [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_run

FUNCTION

  set all options and run a simulation

COPYRIGHT

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

INPUTS

  this=ctqmc
  opt_order=maximal perturbation order to scope
  opt_movie=draw a movie of the simulation
  opt_analysis=compute correlation functions
  opt_check=check fast calculations
  opt_noise=compute noise for green function
  opt_spectra=fourier transform of the time evolution of the number of electrons
  opt_gMove=steps without global move

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

1305 SUBROUTINE Ctqmc_run(this,opt_order,opt_histo,opt_movie,opt_analysis,opt_check,opt_noise,opt_spectra,opt_gMove)
1306 
1307 
1308 #ifdef HAVE_MPI1
1309 include 'mpif.h'
1310 #endif
1311 !Arguments ------------------------------------
1312   TYPE(Ctqmc), INTENT(INOUT)           :: this
1313   INTEGER, OPTIONAL, INTENT(IN   )  :: opt_order
1314   INTEGER, OPTIONAL, INTENT(IN   )  :: opt_histo
1315   INTEGER, OPTIONAL, INTENT(IN   )  :: opt_movie
1316   INTEGER, OPTIONAL, INTENT(IN   )  :: opt_analysis
1317   INTEGER, OPTIONAL, INTENT(IN   )  :: opt_check
1318   INTEGER, OPTIONAL, INTENT(IN   )  :: opt_noise
1319   INTEGER, OPTIONAL, INTENT(IN   )  :: opt_spectra
1320   INTEGER, OPTIONAL, INTENT(IN   )  :: opt_gMove
1321 !Local variables ------------------------------
1322 #ifdef HAVE_MPI
1323   INTEGER                            :: ierr
1324 #endif
1325 !#ifdef CTCtqmc_MOVIE
1326   INTEGER                            :: ilatex
1327   CHARACTER(LEN=4)                   :: Cchar
1328 !#endif
1329   DOUBLE PRECISION                   :: estimatedTime
1330 
1331   IF ( .NOT. this%set  ) &
1332     CALL ERROR("Ctqmc_run : QMC not set up                          ")
1333   IF ( .NOT. this%setU ) &
1334     CALL ERROR("Ctqmc_run : QMC does not have a U this            ")
1335 
1336 
1337 ! OPTIONS of the run
1338   IF ( PRESENT( opt_check ) ) THEN
1339     this%opt_check = opt_check
1340     CALL ImpurityOperator_doCheck(this%Impurity,opt_check)
1341     CALL BathOperator_doCheck(this%Bath,opt_check)
1342   END IF
1343   IF ( PRESENT( opt_movie ) ) &
1344     this%opt_movie = opt_movie
1345   IF ( PRESENT( opt_analysis ) ) &
1346     this%opt_analysis = opt_analysis
1347   IF ( PRESENT ( opt_order ) ) &
1348     this%opt_order = opt_order 
1349   IF ( PRESENT ( opt_histo ) ) &
1350     this%opt_histo = opt_histo 
1351   IF ( PRESENT ( opt_noise ) ) THEN
1352     this%opt_noise = opt_noise 
1353   END IF
1354   IF ( PRESENT ( opt_spectra ) ) &
1355     this%opt_spectra = opt_spectra
1356 
1357   this%modGlobalMove(1) = this%sweeps+1 ! No Global Move
1358   this%modGlobalMove(2) = 0
1359   IF ( PRESENT ( opt_gMove ) ) THEN
1360     IF ( opt_gMove .LE. 0 .OR. opt_gMove .GT. this%sweeps ) THEN
1361       this%modGlobalMove(1) = this%sweeps+1
1362       CALL WARNALL("Ctqmc_run : global moves option is <= 0 or > sweeps/cpu -> No global Moves")
1363     ELSE 
1364       this%modGlobalMove(1) = opt_gMove 
1365     END IF
1366   END IF
1367 
1368   CALL Ctqmc_allocateOpt(this)
1369   
1370 !#ifdef CTCtqmc_MOVIE  
1371   ilatex = 0
1372   IF ( this%opt_movie .EQ. 1 ) THEN
1373     Cchar ="0000"
1374     WRITE(Cchar,'(I4)') this%rank 
1375     ilatex = 87+this%rank
1376     OPEN(UNIT=ilatex, FILE="Movie_"//TRIM(ADJUSTL(Cchar))//".tex")
1377     WRITE(ilatex,'(A)') "\documentclass{beamer}"
1378     WRITE(ilatex,'(A)') "\usepackage{color}"
1379     WRITE(ilatex,'(A)') "\setbeamersize{sidebar width left=0pt}"
1380     WRITE(ilatex,'(A)') "\setbeamersize{sidebar width right=0pt}"
1381     WRITE(ilatex,'(A)') "\setbeamersize{text width left=0pt}"
1382     WRITE(ilatex,'(A)') "\setbeamersize{text width right=0pt}"
1383     WRITE(ilatex,*) 
1384     WRITE(ilatex,'(A)') "\begin{document}"
1385     WRITE(ilatex,*) 
1386   END IF
1387 !#endif
1388 
1389   IF ( this%rank .EQ. 0 ) THEN
1390     WRITE(this%ostream,'(A29)') "Starting QMC (Thermalization)"
1391   END IF
1392   
1393   !=================================
1394   ! STARTING THERMALIZATION 
1395   !=================================
1396   CALL Ctqmc_loop(this,this%thermalization,ilatex)
1397   !=================================
1398   ! ENDING   THERMALIZATION 
1399   !=================================
1400 
1401   estimatedTime = this%runTime
1402 #ifdef HAVE_MPI
1403   CALL MPI_REDUCE(this%runTime, estimatedTime, 1, MPI_DOUBLE_PRECISION, MPI_MAX, &
1404              0, this%MY_COMM, ierr)
1405 #endif
1406 
1407   IF ( this%rank .EQ. 0 ) THEN
1408     WRITE(this%ostream,'(A26,I6,A11)') "Thermalization done in    ", CEILING(estimatedTime), "    seconds"
1409     WRITE(this%ostream,'(A25,I7,A15,I5,A5)') "The QMC should run in    ", &
1410            CEILING(estimatedTime*DBLE(this%sweeps)/DBLE(this%thermalization)),&
1411                         "    seconds on ", this%size, " CPUs"
1412   END IF
1413 
1414   !=================================
1415   ! CLEANING CTQMC          
1416   !=================================
1417   CALL Ctqmc_clear(this)
1418 
1419   !=================================
1420   ! STARTING CTQMC          
1421   !=================================
1422   CALL Ctqmc_loop(this,this%sweeps,ilatex)
1423   !=================================
1424   ! ENDING   CTQMC          
1425   !=================================
1426 
1427   IF ( this%opt_movie .EQ. 1 ) THEN
1428     WRITE(ilatex,*) ""
1429     WRITE(ilatex,'(A14)') "\end{document}"
1430     CLOSE(ilatex)
1431   END IF
1432 
1433   this%done     = .TRUE.
1434 
1435 END SUBROUTINE Ctqmc_run

ABINIT/m_Ctqmc/Ctqmc_setG0wTab [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_setG0wTab

FUNCTION

  Set Gow from input array

COPYRIGHT

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

INPUTS

  this=ctqmc
  Gomega=G0w
  opt_fk=F is already inversed with out iwn

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

814 SUBROUTINE Ctqmc_setG0wTab(this,Gomega,opt_fk)
815 
816 !Arguments ------------------------------------
817   TYPE(Ctqmc), INTENT(INOUT)                      :: this
818   COMPLEX(KIND=8), DIMENSION(:,:), INTENT(IN ) :: Gomega
819   INTEGER                         , INTENT(IN ) :: opt_fk
820 !Local variable -------------------------------
821   DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: F
822 
823   IF ( .NOT. this%para ) &
824     CALL ERROR("Ctqmc_setG0wTab : Ctqmc_setParameters never called    ") 
825 
826   MALLOC(F,(1:this%samples+1,1:this%flavors))
827   CALL Ctqmc_computeF(this,Gomega, F, opt_fk)  ! mu is changed
828   CALL BathOperator_setF(this%Bath, F)
829   !CALL BathOperator_printF(this%Bath)
830   FREE(F)
831 
832   IF ( this%opt_levels .NE. 1 ) THEN ! We compute the mu by hand in computeF
833     CALL ImpurityOperator_setMu(this%Impurity,this%mu)
834   END IF
835 
836   this%inF = .TRUE.
837   this%set = .TRUE. 
838 
839 END SUBROUTINE Ctqmc_setG0wTab

ABINIT/m_Ctqmc/Ctqmc_setMu [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_setMu

FUNCTION

  impose energy levels

COPYRIGHT

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

INPUTS

  this=ctqmc
  levels=energy levels vector

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

SOURCE

1075 SUBROUTINE Ctqmc_setMu(this, levels)
1076 
1077 !Arguments ------------------------------------
1078   TYPE(Ctqmc)                   , INTENT(INOUT) :: this
1079   DOUBLE PRECISION, DIMENSION(:), INTENT(IN   ) :: levels
1080 
1081   IF ( this%flavors .NE. SIZE(levels,1) ) &
1082     CALL WARNALL("Ctqmc_setMu : Taking energy levels from weiss G(iw)")
1083 
1084   this%mu(:)=-levels  ! levels = \epsilon_j - \mu
1085   !this%mu =\tilde{\mu} = \mu -\epsilon_j
1086   CALL ImpurityOperator_setMu(this%Impurity,this%mu)
1087   this%opt_levels = 1
1088 END SUBROUTINE Ctqmc_setMu

ABINIT/m_Ctqmc/Ctqmc_setParameters [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_setParameters

FUNCTION

  set all parameters and operators

COPYRIGHT

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

INPUTS

  this=ctqmc
  buffer=input parameters

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

456 SUBROUTINE Ctqmc_setParameters(this,buffer)
457 
458 !Arguments ------------------------------------
459   TYPE(Ctqmc), INTENT(INOUT)                         :: this
460   DOUBLE PRECISION, DIMENSION(1:9), INTENT(IN   ) :: buffer
461 
462 
463   this%thermalization = INT(buffer(3)) !this%thermalization
464   CALL Ctqmc_setSeed(this,INT(buffer(1)))
465   CALL Ctqmc_setSweeps(this,buffer(2))
466 
467   this%measurements   = INT(buffer(4)) !this%measurements
468   this%flavors        = INT(buffer(5)) !this%flavors
469   this%samples        = INT(buffer(6)) !this%samples
470   this%beta           = buffer(7)      !this%beta
471   this%U              = buffer(8)      !U
472 !  this%mu             = buffer(9)      !this%mu
473   !this%Wmax           = INT(buffer(9)) !Freq
474 !#ifdef CTCtqmc_ANALYSIS
475 !  this%order          = INT(buffer(10)) ! order
476   this%inv_dt         = this%samples / this%beta
477 !#endif
478 
479   !CALL ImpurityOperator_init(this%Impurity,this%flavors,this%beta, this%samples)
480   CALL ImpurityOperator_init(this%Impurity,this%flavors,this%beta)
481   IF ( this%U .GE. 0.d0 ) THEN
482     CALL ImpurityOperator_computeU(this%Impurity,this%U,0.d0)
483     this%setU = .TRUE.
484   END IF
485 !  this%mu = this%mu + this%Impurity%shift_mu
486 
487   CALL BathOperator_init(this%Bath, this%flavors, this%samples, this%beta, INT(buffer(9)))
488 
489   this%para = .TRUE.
490 
491 END SUBROUTINE Ctqmc_setParameters

ABINIT/m_Ctqmc/Ctqmc_setSeed [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_setSeed

FUNCTION

  initialize random number generator

COPYRIGHT

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

INPUTS

  this=ctqmc
  iseed=seed from imput

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

574 SUBROUTINE Ctqmc_setSeed(this,iseed)
575 
576 !Arguments ------------------------------------
577   TYPE(Ctqmc), INTENT(INOUT)           :: this
578   INTEGER  , INTENT(IN   )           :: iseed
579 !Local variables ------------------------------
580   !INTEGER                            :: n
581   !INTEGER                            :: i
582   !INTEGER, DIMENSION(:), ALLOCATABLE :: seed
583 
584 
585   !CALL RANDOM_SEED(size = n)
586   !MALLOC(seed,(n))
587   !seed =  iseed + (/ (i - 1, i = 1, n) /)
588 
589   !CALL RANDOM_SEED(PUT = seed+this%rank)
590 
591   !FREE(seed)
592 
593   this%seed=INT(iseed+this%rank,8)
594 
595 END SUBROUTINE Ctqmc_setSeed

ABINIT/m_Ctqmc/Ctqmc_setSweeps [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_setSweeps

FUNCTION

  set the number of sweeps

COPYRIGHT

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

INPUTS

  this=ctqmc
  sweeps=asked sweeps

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

519 SUBROUTINE Ctqmc_setSweeps(this,sweeps)
520 
521 !Arguments ------------------------------------
522   TYPE(Ctqmc)         , INTENT(INOUT) :: this
523   DOUBLE PRECISION  , INTENT(IN   ) :: sweeps
524 
525   this%sweeps = NINT(sweeps / DBLE(this%size))
526 !  write(6,*) this%sweeps,NINT(sweeps / DBLE(this%size)),ANINT(sweeps/DBLE(this%size))
527   IF ( DBLE(this%sweeps) .NE. ANINT(sweeps/DBLE(this%size)) ) &
528     CALL ERROR("Ctqmc_setSweeps : sweeps is negative or too big     ")
529   IF ( this%sweeps .LT. 2*CTQMC_SLICE1 ) THEN  !202
530     CALL WARNALL("Ctqmc_setSweeps : # sweeps automtically changed     ")
531     this%sweeps = 2*CTQMC_SLICE1
532 !  ELSE IF ( this%sweeps .LT. this%thermalization ) THEN
533 !    CALL WARNALL("Ctqmc_setSweeps : Thermalization > sweeps / cpu -> auto fix")
534 !    this%sweeps = this%thermalization
535   END IF
536   IF ( DBLE(NINT(DBLE(this%sweeps)*DBLE(this%size)/DBLE(CTQMC_SLICE1))) .NE.  &
537   ANINT(DBLE(this%sweeps)*DBLE(this%size)/DBLE(CTQMC_SLICE1)) ) THEN
538     this%modNoise1 = this%sweeps
539   ELSE
540     this%modNoise1    = MIN(this%sweeps,INT(DBLE(this%sweeps)*DBLE(this%size) / DBLE(CTQMC_SLICE1))) !101
541   END IF
542   this%modNoise2    = MAX(this%modNoise1 / CTQMC_SLICE2, 1)   ! 100
543 !  this%modGlobalMove(1) = this%thermalization / 10 + 1
544 !  this%modGlobalMove(2) = 0
545 
546 END SUBROUTINE Ctqmc_setSweeps

ABINIT/m_Ctqmc/Ctqmc_setU [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_setU

FUNCTION

  set the interaction this

COPYRIGHT

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

INPUTS

  this=ctqmc
  matU=interaction this

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

898 SUBROUTINE Ctqmc_setU(this,matU)
899 
900 !Arguments ------------------------------------
901   TYPE(Ctqmc), INTENT(INOUT) :: this
902 !Local variables ------------------------------
903   DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN) :: matU
904 
905   IF ( SIZE(matU) .NE. this%flavors*this%flavors ) &
906     CALL ERROR("Ctqmc_setU : Wrong interaction this (size)        ")
907 
908   CALL ImpurityOperator_setUmat(this%Impurity, matU)
909   this%setU = .TRUE.
910 END SUBROUTINE Ctqmc_setU

ABINIT/m_Ctqmc/Ctqmc_symmetrizeGreen [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_symmetrizeGreen

FUNCTION

  optionnaly symmetrize the green functions

COPYRIGHT

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

INPUTS

  this=ctqmc
  syms=weight factors

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

SOURCE

2564 SUBROUTINE Ctqmc_symmetrizeGreen(this, syms)
2565 
2566 !Arguments ------------------------------------
2567   TYPE(Ctqmc)                     , INTENT(INOUT) :: this
2568   DOUBLE PRECISION, DIMENSION(:,:), INTENT(IN   ) :: syms
2569 !Local variables ------------------------------
2570   INTEGER :: iflavor1
2571   INTEGER :: iflavor2
2572   INTEGER :: flavors
2573   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: green_tmp
2574   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:  ) :: n_tmp
2575 
2576   flavors = this%flavors
2577   IF ( SIZE(syms,1) .NE. flavors .OR. SIZE(syms,2) .NE. flavors ) THEN
2578     CALL WARNALL("Ctqmc_symmetrizeGreen : wrong opt_sym -> not symmetrizing")
2579     RETURN
2580   END IF
2581  
2582   MALLOC(green_tmp,(1:this%samples+1,flavors))
2583   green_tmp(:,:) = 0.d0
2584   MALLOC(n_tmp,(1:flavors))
2585   n_tmp(:) = 0.d0
2586   DO iflavor1=1, flavors
2587     DO iflavor2=1,flavors
2588       green_tmp(:,iflavor1) = green_tmp(:,iflavor1) &
2589                              + syms(iflavor2,iflavor1) * this%Greens(iflavor2)%oper(:)
2590       n_tmp(iflavor1) = n_tmp(iflavor1) &
2591                              + syms(iflavor2,iflavor1) * this%measN(1,iflavor2)
2592     END DO
2593   END DO
2594   DO iflavor1=1, flavors
2595     this%Greens(iflavor1)%oper(:) = green_tmp(:,iflavor1)
2596     this%measN(1,iflavor1)          = n_tmp(iflavor1)
2597   END DO
2598   FREE(green_tmp)
2599   FREE(n_tmp)
2600 END SUBROUTINE Ctqmc_symmetrizeGreen

ABINIT/m_Ctqmc/Ctqmc_tryAddRemove [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_tryAddRemove

FUNCTION

  Try to add or remove a segment and an anti-segment

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

  updated=something changed

SIDE EFFECTS

NOTES

SOURCE

1700 SUBROUTINE Ctqmc_tryAddRemove(this,updated)
1701 
1702 !Arguments ------------------------------------
1703   TYPE(Ctqmc)             , INTENT(INOUT) :: this
1704 !  TYPE(BathOperator)    , INTENT(INOUT) :: Bath 
1705 !  TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity 
1706   LOGICAL               , INTENT(  OUT) :: updated
1707 !Local variables ------------------------------
1708   INTEGER                               :: position
1709   INTEGER         , DIMENSION(1:2)     :: nature ! -2 for antiseg and 1 for seg
1710   INTEGER                               :: i! -2 for antiseg and 1 for seg
1711   DOUBLE PRECISION                      :: action
1712   DOUBLE PRECISION                      :: beta
1713   DOUBLE PRECISION                      :: time1
1714   DOUBLE PRECISION                      :: time2
1715   DOUBLE PRECISION                      :: time_avail
1716   DOUBLE PRECISION                      :: det_ratio
1717   DOUBLE PRECISION                      :: imp_trace
1718   DOUBLE PRECISION                      :: signe
1719   DOUBLE PRECISION                      :: tail
1720   DOUBLE PRECISION, DIMENSION(1:2)      :: CdagC_1
1721 
1722   IF ( .NOT. this%set ) &
1723     CALL ERROR("Ctqmc_trySegment : QMC not set                       ")
1724 
1725   nature(1) = CTQMC_SEGME
1726   nature(2) = CTQMC_ANTIS
1727   beta      = this%beta
1728 
1729   updated = .FALSE.
1730   tail  = DBLE(this%Impurity%particles(this%Impurity%activeFlavor)%tail)
1731 
1732 
1733   DO i = 1, 2
1734     signe = SIGN(1.d0,DBLE(nature(i))) 
1735 
1736     !CALL RANDOM_NUMBER(action)
1737     CALL OurRng(this%seed,action)
1738 
1739     IF ( action .LT. .5d0 ) THEN ! Ajout de segment
1740       !CALL RANDOM_NUMBER(time1)
1741       CALL OurRng(this%seed,time1)
1742       time1 = time1 * beta
1743       time_avail = ImpurityOperator_getAvailableTime(this%Impurity,time1,position) * signe
1744       IF ( time_avail .GT. 0.d0 ) THEN
1745         !CALL RANDOM_NUMBER(time2)
1746         CALL OurRng(this%seed,time2)
1747         IF ( time2 .EQ. 0.d0 ) THEN
1748           CALL OurRng(this%seed,time2) ! Prevent null segment
1749         END IF
1750         time2     = time1 + time2 * time_avail
1751         CdagC_1(Cdag_) = ((1.d0+signe)*time1+(1.d0-signe)*time2)*0.5d0
1752         CdagC_1(C_   ) = ((1.d0+signe)*time2+(1.d0-signe)*time1)*0.5d0
1753         det_ratio = BathOperator_getDetAdd(this%Bath,CdagC_1,position,this%Impurity%particles(this%Impurity%activeFlavor))
1754         imp_trace = ImpurityOperator_getTraceAdd(this%Impurity,CdagC_1)
1755         !CALL RANDOM_NUMBER(time1)
1756         CALL OurRng(this%seed,time1)
1757         IF ( det_ratio*imp_trace .LT. 0.d0 ) THEN
1758           this%stats(nature(i)+CTQMC_DETSI) = this%stats(nature(i)+CTQMC_DETSI) + 1.d0
1759         END IF
1760         IF ( (time1 * (tail + 1.d0 )) &
1761              .LT. (beta * time_avail * det_ratio * imp_trace ) ) THEN
1762           CALL ImpurityOperator_add(this%Impurity,CdagC_1,position)
1763           CALL BathOperator_setMAdd(this%bath,this%Impurity%particles(this%Impurity%activeFlavor))
1764           this%stats(nature(i)+CTQMC_ADDED) = this%stats(nature(i)+CTQMC_ADDED)  + 1.d0
1765           updated = .TRUE. .OR. updated
1766           tail = tail + 1.d0
1767         END IF 
1768       END IF 
1769 
1770     ELSE ! Supprimer un segment
1771       IF ( tail .GT. 0.d0 ) THEN
1772         !CALL RANDOM_NUMBER(time1)
1773         CALL OurRng(this%seed,time1)
1774         position = INT(((time1 * tail) + 1.d0) * signe )
1775         time_avail = ImpurityOperator_getAvailedTime(this%Impurity,position)
1776         det_ratio  = BathOperator_getDetRemove(this%Bath,position)
1777         imp_trace  = ImpurityOperator_getTraceRemove(this%Impurity,position)
1778         !CALL RANDOM_NUMBER(time1)
1779         CALL OurRng(this%seed,time1)
1780         IF ( det_ratio * imp_trace .LT. 0.d0 ) THEN
1781           this%stats(nature(i)+CTQMC_DETSI) = this%stats(nature(i)+CTQMC_DETSI) + 1.d0
1782         END IF
1783         IF ( (time1 * beta * time_avail ) &
1784              .LT. (tail * det_ratio * imp_trace) ) THEN
1785           CALL ImpurityOperator_remove(this%Impurity,position)
1786           CALL BathOperator_setMRemove(this%Bath,this%Impurity%particles(this%Impurity%activeFlavor))
1787           this%stats(nature(i)+CTQMC_REMOV) = this%stats(nature(i)+CTQMC_REMOV)  + 1.d0
1788           updated = .TRUE. .OR. updated
1789           tail = tail -1.d0
1790         END IF
1791       END IF
1792     END IF
1793   END DO
1794 END SUBROUTINE Ctqmc_tryAddRemove

ABINIT/m_Ctqmc/Ctqmc_trySwap [ Functions ]

[ Top ] [ Functions ]

NAME

  Ctqmc_trySwap

FUNCTION

  try a global move (swap to flavors)

COPYRIGHT

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

INPUTS

  this=ctqmc

OUTPUT

  flav_i=first flavor swaped
  flav_j=second flavor swaped

SIDE EFFECTS

NOTES

SOURCE

1823 SUBROUTINE Ctqmc_trySwap(this,flav_i,flav_j)
1824 
1825 !Arguments ------------------------------------
1826   TYPE(Ctqmc)           , INTENT(INOUT) :: this
1827 !  TYPE(BathOperator)    , INTENT(INOUT) :: Bath 
1828 !  TYPE(ImpurityOperator), INTENT(INOUT) :: Impurity 
1829   INTEGER               , INTENT(  OUT) :: flav_i
1830   INTEGER               , INTENT(  OUT) :: flav_j
1831 !Local variables ------------------------------
1832   INTEGER :: flavor_i
1833   INTEGER :: flavor_j
1834   DOUBLE PRECISION :: rnd
1835   DOUBLE PRECISION :: lengthi
1836   DOUBLE PRECISION :: lengthj
1837   DOUBLE PRECISION :: overlapic1
1838   DOUBLE PRECISION :: overlapjc1
1839   DOUBLE PRECISION :: overlapic2
1840   DOUBLE PRECISION :: overlapjc2
1841   DOUBLE PRECISION :: detic1
1842   DOUBLE PRECISION :: detjc1
1843   DOUBLE PRECISION :: detic2
1844   DOUBLE PRECISION :: detjc2
1845   DOUBLE PRECISION :: det_ratio
1846   DOUBLE PRECISION :: local_ratio
1847 
1848   !CALL RANDOM_NUMBER(rnd)
1849   CALL OurRng(this%seed,rnd)
1850   flavor_i = NINT(rnd*DBLE(this%flavors-1.d0))+1
1851   !CALL RANDOM_NUMBER(rnd)
1852   CALL OurRng(this%seed,rnd)
1853   flavor_j = NINT(rnd*DBLE(this%flavors-1.d0))+1
1854   
1855   flav_i = 0
1856   flav_j = 0
1857 
1858   IF ( flavor_i .NE. flavor_j ) THEN
1859     ! On tente d'intervertir i et j
1860     ! Configuration actuelle :
1861     this%modGlobalMove(2) = this%modGlobalMove(2)+1
1862     detic1     = BathOperator_getDetF(this%Bath,flavor_i)
1863     detjc1     = BathOperator_getDetF(this%Bath,flavor_j)
1864     lengthi    = ImpurityOperator_measN(this%Impurity,flavor_i)
1865     lengthj    = ImpurityOperator_measN(this%Impurity,flavor_j)
1866     overlapic1 = ImpurityOperator_overlapFlavor(this%Impurity,flavor_i)
1867     overlapjc1 = ImpurityOperator_overlapFlavor(this%Impurity,flavor_j)
1868     ! Configuration nouvelle :
1869     detic2     = BathOperator_getDetF(this%Bath,flavor_i,this%Impurity%particles(flavor_j))
1870     detjc2     = BathOperator_getDetF(this%Bath,flavor_j,this%Impurity%particles(flavor_i))
1871     ! lengths unchanged
1872     overlapic2 = ImpurityOperator_overlapSwap(this%Impurity,flavor_i,flavor_j)
1873     overlapjc2 = ImpurityOperator_overlapSwap(this%Impurity,flavor_j,flavor_i)
1874 
1875 !    IF ( detic1*detjc1 .EQ. detic2*detjc2 ) THEN
1876 !      det_ratio = 1.d0
1877 !    ELSE IF ( detic1*detjc1 .EQ. 0.d0 ) THEN
1878 !      det_ratio = detic2*detjc2 ! evite de diviser par 0 si pas de segment
1879 !    ELSE
1880       det_ratio = detic2*detjc2/(detic1*detjc1)
1881 !    END IF
1882     local_ratio = DEXP(-overlapic2*overlapjc2+overlapic1*overlapjc1 &
1883                       +(lengthj-lengthi)*(this%mu(flavor_i)-this%mu(flavor_j)))
1884     ! Wloc = exp(muN-Uo)
1885     !CALL RANDOM_NUMBER(rnd)
1886     CALL OurRng(this%seed,rnd)
1887     IF ( rnd .LT. local_ratio*det_ratio ) THEN ! swap accepted
1888       CALL ImpurityOperator_swap(this%Impurity, flavor_i,flavor_j)
1889       CALL BathOperator_swap    (this%Bath    , flavor_i,flavor_j)
1890       this%swap = this%swap + 1.d0
1891       flav_i = flavor_i
1892       flav_j = flavor_j
1893 !    ELSE
1894 !      CALL WARN("Swap refused")
1895 !      WRITE(this%ostream,'(6E24.14)') local_ratio, det_ratio, detic1, detjc1, detic2, detjc2
1896     END IF
1897   END IF
1898 
1899 END SUBROUTINE Ctqmc_trySwap

m_Ctqmc/Ctqmc [ Types ]

[ Top ] [ m_Ctqmc ] [ Types ]

NAME

  Ctqmc

FUNCTION

  This structured datatype contains the necessary data

COPYRIGHT

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

SOURCE

 73 TYPE, PUBLIC :: Ctqmc
 74 
 75   LOGICAL _PRIVATE :: init = .FALSE.
 76 ! Flag: is MC initialized
 77 
 78   LOGICAL _PRIVATE :: set  = .FALSE.
 79 ! Flag: ??
 80 
 81   LOGICAL _PRIVATE :: setU = .FALSE.
 82 ! Flag: is U Set ?
 83 
 84   LOGICAL _PRIVATE :: inF  = .FALSE.
 85 ! Flag: is hybridization fct in input ?
 86 
 87   LOGICAL _PRIVATE :: done = .FALSE.
 88 ! Flag: is MC terminated ?
 89 
 90   LOGICAL _PRIVATE :: para = .FALSE.
 91 ! Flag:  do we have parameters in input
 92 
 93   LOGICAL _PRIVATE :: have_MPI = .FALSE.
 94 ! Flag: 
 95 
 96   INTEGER _PRIVATE :: opt_movie = 0
 97 !
 98 
 99   INTEGER _PRIVATE :: opt_analysis = 0
100 ! correlations 
101 
102   INTEGER _PRIVATE :: opt_check = 0
103 ! various check 0
104 ! various check 1 impurity
105 ! various check 2 bath
106 ! various check 3 both
107 
108   INTEGER _PRIVATE :: opt_order = 0
109 ! nb of segments max for analysis
110 
111   INTEGER _PRIVATE :: opt_histo = 0
112 ! Enable histograms
113 
114   INTEGER _PRIVATE :: opt_noise = 0
115 ! compute noise
116 
117   INTEGER _PRIVATE :: opt_spectra = 0
118 ! markov chain FT (correlation time)
119 
120   INTEGER _PRIVATE :: opt_levels = 0
121 ! do we have energy levels
122 
123   INTEGER _PRIVATE :: flavors
124 !
125 
126   INTEGER _PRIVATE :: measurements
127 ! nb of measure in the MC
128 
129   INTEGER _PRIVATE :: samples
130 ! nb of L points
131 
132   INTEGER(8) _PRIVATE :: seed
133 !
134 
135   INTEGER _PRIVATE :: sweeps
136 !
137 
138   INTEGER _PRIVATE :: thermalization
139 !
140 
141   INTEGER _PRIVATE :: ostream
142 ! output file
143 
144   INTEGER _PRIVATE :: istream
145 ! input file
146 
147   INTEGER _PRIVATE :: modNoise1
148 ! measure the noise each modNoise1
149 
150   INTEGER _PRIVATE :: modNoise2
151 ! measure the noise each modNoise2
152 
153   INTEGER _PRIVATE :: activeFlavor
154 ! orbital on which one do sth now
155 
156   INTEGER, DIMENSION(1:2) _PRIVATE :: modGlobalMove
157 ! 1: gloabl move each modglobalmove(1)
158 ! 2: we have done modglobalmove(2) for two different orbitals.
159 
160   INTEGER _PRIVATE :: Wmax
161 ! Max freq for FT
162 
163   DOUBLE PRECISION, DIMENSION(1:6) _PRIVATE :: stats
164 ! to now how many negative determinant, antisegments,seeme.e.twfs...j
165 
166   DOUBLE PRECISION _PRIVATE :: swap
167 ! nb of successfull GM
168 
169   INTEGER _PRIVATE :: MY_COMM
170 ! 
171 
172   INTEGER _PRIVATE :: rank
173 !
174 
175   INTEGER _PRIVATE :: size
176 ! size of MY_COMM
177 
178   DOUBLE PRECISION _PRIVATE :: runTime ! time for the run routine
179 !  
180 
181   DOUBLE PRECISION _PRIVATE :: beta
182 !
183 
184   DOUBLE PRECISION _PRIVATE :: U
185 
186   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) _PRIVATE :: mu
187 ! levels
188 
189   TYPE(GreenHyb)  , ALLOCATABLE, DIMENSION(:    ) _PRIVATE :: Greens 
190 
191   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:  ) _PRIVATE :: measN 
192 ! measure of occupations (3or4,flavor) 
193 
194   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:  ) _PRIVATE :: measDE
195 !  (flavor,flavor) double occupancies
196 !  (1,1): total energy of correlation.
197 
198   DOUBLE PRECISION _PRIVATE :: a_Noise
199 ! Noise a exp (-bx) for the  noise
200 
201   DOUBLE PRECISION _PRIVATE :: b_Noise
202 ! Noise a exp (-bx) for the  noise
203 
204   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: abNoiseG   !(ab,tau,flavor)
205 ! Noise but for G
206 
207   TYPE(Vector)             , DIMENSION(1:2) _PRIVATE :: measNoise 
208   TYPE(Vector), ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: measNoiseG       !(tau,flavor,mod) 
209 ! accumulate each value relataed to measurenoise 1 2
210 
211   DOUBLE PRECISION _PRIVATE                            :: inv_dt
212 ! 1/(beta/L)
213 
214   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:  ) _PRIVATE :: measPerturbation 
215 ! opt_order,nflavor
216 
217   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: occup_histo_time
218 ! nflavor
219 
220   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) _PRIVATE :: measCorrelation 
221 ! segment,antisegment,nflavor,nflavor
222 
223   DOUBLE PRECISION _PRIVATE :: errorImpurity
224 ! check 
225 
226   DOUBLE PRECISION _PRIVATE :: errorBath
227 ! for check
228 
229   TYPE(BathOperator) _PRIVATE              :: Bath
230 
231   TYPE(ImpurityOperator) _PRIVATE          :: Impurity
232 
233   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) _PRIVATE :: density
234 
235 END TYPE Ctqmc