core dump in Trace.F90 due to subscript greater than upper bound

[yuechm]

Dear developer,

Recently I installed the w2dynamics and did a worm measurement test of Green's function in 1 band Bethe lattice. It seems that there is a bug in Trace.F90 in the subroutine get_LogDetRatPairAdd. There is a segment fault (core dump) reported after the running. To find what's going on, I recompiled the code with the debug flag "-warn all, -check all, -traceback". The error is as following

forrtl: severe (408): fort: (2): Subscript #1 of the array R has value 2 which is greater than the upper bound of 1

Image PC Routine Line Source
libifcoremt.so.5 00007F330C818143 for_emit_diagnost Unknown Unknown
CTQMC.so 00007F330F93767A mtrace_mp_get_log 2784 Trace.F90
CTQMC.so 00007F330F68F318 Unknown Unknown Unknown
CTQMC.so 00007F330F6AE625 Unknown Unknown Unknown

My Trace.F90 modified to is
2769 do while (associated(Element))
2770 if(Element%has_hyb) then
2771 if(Element%tau.ne.ElementS%tau)then
2772 print *, "Element%CA ElementS%CA", Element%CA, ElementS%CA
2773 print *, "Element%CA.ne.ElementS%CA", Element%CA.ne.ElementS%CA
2774
2775 print *, "Element%CA ElementE%CA", Element%CA, ElementE%CA
2776 print *, "Element%CA.ne.ElementE%CA", Element%CA.ne.ElementE%CA
2777
2778 if(Element%CA.ne.ElementS%CA.and..not.associated(Element,ElementE))then
2779 print *, "now iQ", iQ
2780 Q(iQ)=get_HybrF(this,Element%tau-ElementS%tau,Element%Orbital,Element%Spin)
2781 iQ=iQ+1
2782 elseif(Element%CA.ne.ElementE%CA.and..not.associated(Element,ElementS))then
2783 print *, "now iR", iR
2784 R(iR)=get_HybrF(this,ElementE%tau-Element%tau,Element%Orbital,Element%Spin)
2785 iR=iR+1
2786 endif
2787 endif
2788 endif
2789 Element=>Element%fnext
2790 enddo

The last output is
size Q 1
size R 1
Element%CA ElementS%CA 1 1
Element%CA.ne.ElementS%CA F
Element%CA ElementE%CA 1 2
Element%CA.ne.ElementE%CA T
now iR 1
Element%CA ElementS%CA 2 1
Element%CA.ne.ElementS%CA T
Element%CA ElementE%CA 2 2
Element%CA.ne.ElementE%CA F
now iQ 1
Element%CA ElementS%CA 1 1
Element%CA.ne.ElementS%CA F
Element%CA ElementE%CA 1 2
Element%CA.ne.ElementE%CA T
now iR 2

forrtl: severe (408): fort: (2): Subscript #1 of the array R has value 3 which is greater than the upper bound of 2

Image PC Routine Line Source
libifcoremt.so.5 00007EFE53257143 for_emit_diagnost Unknown Unknown
CTQMC.so 00007EFE5637667A mtrace_mp_get_log 2784 Trace.F90
CTQMC.so 00007EFE560CE318 Unknown Unknown Unknown
CTQMC.so 00007EFE560C9AA5 Unknown Unknown Unknown
CTQMC.so 00007EFE560C7DFE Unknown Unknown Unknown
CTQMC.so 00007EFE560B2FF8 Unknown Unknown Unknown
libpython2.7.so.1 00007EFE6AFE4B73 PyObject_Call Unknown Unknown
libpython2.7.so.1 00007EFE6B07B119 PyEval_EvalFrameE Unknown Unknown
libpython2.7.so.1 00007EFE6B080A99 PyEval_EvalCodeEx Unknown Unknown
libpython2.7.so.1 00007EFE6B07DF68 PyEval_EvalFrameE Unknown Unknown
libpython2.7.so.1 00007EFE6B080A99 PyEval_EvalCodeEx Unknown Unknown
libpython2.7.so.1 00007EFE6B080CBA PyEval_EvalCode Unknown Unknown
libpython2.7.so.1 00007EFE6B09A01D Unknown Unknown Unknown
libpython2.7.so.1 00007EFE6B09B1C8 PyRun_FileExFlags Unknown Unknown
libpython2.7.so.1 00007EFE6B09C3E8 PyRun_SimpleFileE Unknown Unknown
libpython2.7.so.1 00007EFE6B0AE67C Py_Main Unknown Unknown
libc-2.17.so 00007EFE6A294C05 __libc_start_main Unknown Unknown
python2.7 000055DE155E707F Unknown Unknown Unknown

The NSeed I provide to DMFT.py is
Nseed = 1993161978+mpi_rank
print "AAAATTTTENNTION! SEED IS CONST FOR ANY RANK!!!", Nseed
#Nseed = cfg["QMC"]["NSeed"] + mpi_rank

My Parameters.in is

[General]
DOS = Bethe # parse the hamiltonian defined previously
half-bandwidth = 2., # half bandwidth of Bethe DOS
beta =10 # inverse temperature
NAt = 1 # no. of atoms; identical to [Atoms] subsections
Uw=0
mu = 0.0 # chemical potential
totdens = 1.0 # targeted fillings
EPSN = 0.001 # chemical potential search, 0-> disabled
DMFTsteps = 0 # no. of DMFT steps
StatisticSteps = 0 # stat iter
WormSteps = 1 # no. of worm
magnetism = para # paramagnetic calculation
mixing = 0.5 # mixing of old and new self-energy; 0 means 100% of the new self energy
FileNamePrefix = U_0.0 # output filename, time stamp appended
fileold = asdf # filename of hdf5-file a calculation is continued
readold = 0 # iteration where run is continued; 0 -> new run

[Atoms]
[[1]]
Nd = 1
Hamiltonian = Density # local interaction (Density,Kanamori,Coulomb)
Udd = 2.0 # intra-orbital Coulomb U
Jdd = 0.0 # Hund coupling J
Vdd = 0.0 # inter-orbital Coulomb V
QuantumNumbers = Nt Szt Azt # conserved quantities; Nt = number of electrons,
# Szt = spin in z-direction, Azt = densit, Qzt PS quantum number
[QMC]
NSeed = 43890 # you may delete it
Nwarmups = 100000 # no. of thermalisation steps before measurement
Nmeas = 1000000 # no. of measurement steps
NCorr = 20 # auto-correlation estimate
PercentageGlobalMove = 0.005 # share of global updates
Ntau = 1000 # no. of imaginary-time bins for the 1P GF
Nftau = 2001 # no. of imaginary-time points for the hybridisation function
Niw = 1000 # no. of positive Matsubara frequencies for the 1P GF; note that Nftau > 2*Niw
NLegMax = 30 # no. of Legendre coefficient for the 1P GF
NLegOrder = 30 # no. of Legendre coefficient for the 1P GF used in DMFT loop
ReuseMCConfig = 0

PercentageWormInsert = 0.3 # percentage of worm insert
PercentageWormReplace= 0.1 # percentage of worm remove
WormMeasGtau = 1 # gtau worm
WormEta = 1
WormSearchEta = 0
Nwarmups2Plus = 5e5

Could you please help on checking the "Parameters.in" file for me? Is there something wrong with this input file? If nothing wrong, maybe there is really a bug in Trace/get_LogDetRatPairAdd

Thanks!

Changming Yue
Postdoc at the University of Fribourg.

[alexkowalski]

Dear Changming Yue,

thank you for reporting this bug. If you only try to measure the diagonal components of the Green's function, you should not have a problem. To do that, you can either set WormComponents = 1, 4 or use the automatic eta search by setting WormSearchEta = 1 instead of providing your own eta. This will automatically skip the components that are zero, but I am not sure if you might not still encounter this bug in that case, because after a short look at the relevant code it seems to me that the eta search might result in a call that would cause get_LogDetRatPairAdd to crash as well, although it did not do that when I tried it. You might also want to enable the eta search even if you select the non-zero components explicitly because your eta of 1 would probably cause unbalanced sampling of the Z and worm space and results of worse statistical quality than you could get in the same time with a well chosen eta (which would be about 0.065 in this calculation according to my very short runs). If you want to do different calculations later, you can use the function component2index_general in w2dyn.auxiliaries.compoundIndex to get the right indices that you need to provide in WormComponents to sample the right components.

Best regards,
Alexander Kowalski