This measures connected propagotors and disconnected loops for singlet
mesons with up/down(degenerate) and strange quarks.
program main_singlet
!
! This measures connected propagotors and disconnected loops for singlet mesons
! with up/down(degenerate) and strange quarks.
!
use comlib
use file_tools_class
use lattice_class
use noisy_wavefunc_class
use fftw3_class
use timer_class
use wavefunc_class
use quark_prop_class
use hadron_class
use main_singlet_class
implicit none
character(len=CHARLEN), parameter :: prog_ver="$Id: main_singlet.F90,v 1.7 2011/10/31 07:57:06 ishikawa Exp $"
type(main_singlet_obj) :: main
type(quark_prop), allocatable :: qprop
type(quark_prop), allocatable :: tmpqp
type(src_wavefunc_obj), allocatable, target :: Z2_noise_wf(:)
type(snk_wavefunc_obj), allocatable, target :: snk_wf(:)
type(meson_prop_obj) :: meson
type(meson_prop_obj) :: meson_disc(2)
type(hadrons_obj) :: hadrons
character(len=CHARLEN) :: str(2),char
integer :: traj
integer :: isrc,iqk,iwf,itt_src
integer :: ifm,iout,inoise
integer :: wf_table(2,2)
wf_table(1,1)=1 ! up/down with local wavefunc
wf_table(1,2)=2 ! up/down with smeared wavefunc (param 1)
wf_table(2,1)=1 ! strange with local wavefunc
wf_table(2,2)=3 ! strange with smeared wavefunc (param 2)
call new(main)
if (nodeid == 0) then
write(*,'(80("-"))')
write(*,'(1X,A," for ",A)')TRIM(prog_ver),TRIM(_VERSION_)
write(*,'(80("-"))')
endif
call new(hadrons)
iout= search_free_file_unit()
call read(iout,main)
close(iout)
call print(main)
allocate(qprop)
call new(qprop)
!=================================
! prepare sink wave functions FFT
!=================================
if (nodeid==0) then
write(*,'("")')
write(*,'(80("="))')
write(*,'("Begin FFT Sink wavefunction")')
endif
allocate(snk_wf(3))
do isrc=1,3
call new(snk_wf(isrc)) ! isrc=1 local-smr, =2 up-down-smr, =3 strange-smr
enddo
do iqk =1,2
do isrc=1,2
if (wf_table(iqk,isrc) /= 1) then
snk_wf(wf_table(iqk,isrc))%param = main%hadrons%wavefunc(wf_table(iqk,isrc))%param
call get(snk_wf(wf_table(iqk,isrc)))
call fft_wavefunc(snk_wf(wf_table(iqk,isrc)))
endif
enddo
enddo
if (nodeid==0) then
write(*,'("End FFT Sink wavefunction")')
write(*,'("")')
endif
do traj=main%hadrons%ini_traj,main%hadrons%end_traj,main%hadrons%traj_skip
!================
! read config
!================
call read(traj,main%hadrons%config_path,main%hadrons%config_name, main%hadrons%config_prof,main%hadrons%u)
meson%itraj = traj
meson%path = main%hadrons%hadron_path
meson%config_prof => main%hadrons%config_prof
!!====================================================
!! corrdinate origin shift. (for trans. sym. test)
! call shift(main%hadrons%u, &
! & main%hadrons%wavefunc(1)%itx0, &
! & main%hadrons%wavefunc(1)%ity0, &
! & main%hadrons%wavefunc(1)%itz0)
!!====================================================
!goto 1000 ! skip connected
!====================================================
! Connected propagator
!====================================================
!===============================
! multiple time slice for source
!===============================
do itt_src=1,NTT,NTT/4
! do itt_src=1,1
write(char,'(".",I3.3)')itt_src
meson%fname = "chpropM"//TRIM(ADJUSTL(char)) ! connected meson prop
call new(meson)
if (nodeid==0) then
write(*,'("")')
write(*,'(80("="))')
write(*,'("Begin Connected part, itt_src=",I3)')itt_src
endif
ifm=1
do iqk =1,2
do isrc=1,2 ! isrc=1 for local, 2 for smeared source
call new(qprop)
!===============================
! set mass parameters
!===============================
qprop%parameter => main%hadrons%quark_parameter(iqk)
!===============================
! set wavefunction parameters
!===============================
!------------
! source wf
!------------
qprop%wavefunc_sous => main%hadrons%wavefunc(wf_table(iqk,isrc))
qprop%wavefunc_sous%itt0 = itt_src
qprop%parameter%smear_sous = isrc-1
!------------
! local sink
!------------
qprop%wavefunc_sink => snk_wf(wf_table(iqk,1))
qprop%parameter%smear_sink = 0
if (nodeid==0) then
write(*,'("")')
write(*,'(80("="))')
write(*,'("Begin ==== Quark #",I3," Source #",I3)')iqk,isrc
endif
call print(qprop)
if (nodeid==0) then
write(*,'("")')
endif
!===============================
! generate source wavefunction
!===============================
call get(qprop%wavefunc_sous)
!===================================================
! compute smared source-local sink quark propagator
!===================================================
call set_source(qprop)
call get_prop(qprop,main%hadrons%u)
deallocate(qprop%wavefunc_sous%wavefunc)
if (nodeid==0) then
write(*,'("")')
endif
call print_statistics(qprop)
!==============================================================
! local sink meson propagator
!==============================================================
if (nodeid==0) then
write(*,*)
write(*,'(80("="))')
write(*,'(" Connected Meson propagator local sink")')
write(*,*)
endif
call get_prop(meson,qprop,qprop)
call save(meson,ifm)
ifm=ifm+1
!==============================================================
! double smeared sink meson propagator
!==============================================================
if (nodeid==0) then
write(*,'(80("="))')
write(*,'(" Connected Meson propagator smeared sink")')
write(*,*)
write(*,'(" Smear sink")')
endif
allocate(tmpqp)
call new(tmpqp)
call assign(tmpqp,qprop)
tmpqp%wavefunc_sink => snk_wf(wf_table(iqk,2))
tmpqp%parameter%smear_sink = 1
call fft_qprop(tmpqp%prop)
call convolution(snk_wf(wf_table(iqk,2))%wavefunc,tmpqp%prop)
call invfft_qprop(tmpqp%prop)
call get_prop(meson,tmpqp,tmpqp)
call delete(tmpqp)
deallocate(tmpqp)
call save(meson,ifm)
ifm=ifm+1
call delete(qprop)
if (nodeid==0) then
write(*,'("")')
write(*,'("End ==== Quark #",I3," Source #",I3)')iqk,isrc
endif
enddo
enddo
call delete(meson)
if (nodeid==0) then
write(*,'("End Connected part, itt_src=",I3)')itt_src
endif
enddo
1000 continue
!call delete(main)
!stop
!======================================================================
! Disconnected part
!======================================================================
if (nodeid==0) then
write(*,'("")')
write(*,'(80("="))')
write(*,'("Begin Disconnected part")')
endif
allocate(Z2_noise_wf(2))
str(1) = " Sink loop"
str(2) = " Source loop"
do inoise=1,main%nnoise ! noise loop
call new(Z2_noise_wf(1)) ! for snk side loop
call new(Z2_noise_wf(2)) ! for src side loop
meson_disc(1) = meson ! for snk side loop
meson_disc(2) = meson ! for src side loop
meson_disc(1)%fname = "singlet_snk_dhpropM"
meson_disc(2)%fname = "singlet_src_dhpropM"
call new_meson_dloop(meson_disc(1),inoise)
call new_meson_dloop(meson_disc(2),inoise)
!===============================================
! set Z2 noise source to take the trace prop
!===============================================
call set_noise(Z2_noise_wf(1),m_rand) ! for the sink side trace
call set_noise(Z2_noise_wf(2),m_rand) ! for the source side trace
ifm=1
do iqk =1,2
do isrc=1,2 ! 1 for snk loop, 2 for src loop
if (nodeid==0) then
write(*,'("")')
write(*,'(80("="))')
write(*,'("Begin ==== Quark #",I3," Z2 noise #",I3,A)')iqk,inoise,TRIM(str(isrc))
endif
call new(qprop)
!===============================
! set mass parameters
!===============================
qprop%parameter => main%hadrons%quark_parameter(iqk)
!===============================================
! compute quark propagator with Z2 noise source
!===============================================
qprop%wavefunc_sous => Z2_noise_wf(isrc)
qprop%parameter%smear_sous = 1
qprop%wavefunc_sink => snk_wf(wf_table(iqk,1))
qprop%parameter%smear_sink = 0
call print(qprop)
if (nodeid==0) write(*,'("")')
call set_source(qprop) ! set Z2 noise source
call get_prop(qprop,main%hadrons%u)
if (nodeid==0) write(*,'("")')
call print_statistics(qprop)
!=============================================================
! Compute quark trace loops for disconnected meson propagator
!=============================================================
if (nodeid==0) write(*,'("Start Disconnected Quark Trace Loop")')
!================================
! no smearing
!================================
iwf=0
if (nodeid==0) then
write(*,*)
write(*,'(80("-"))')
write(*,'("Start",A," without smearing")')TRIM(str(isrc))
write(*,*)
endif
call get_meson_dloop(meson_disc(isrc),Z2_noise_wf(isrc),qprop)
call save_meson_dloop(meson_disc(isrc),inoise,ifm,iwf)
ifm=ifm+1
!================================
! doubly smeared
!================================
if (nodeid==0) then
write(*,*)
write(*,'(80("-"))')
write(*,'("Start",A," with double smearing")')TRIM(str(isrc))
write(*,*)
endif
iwf=2
allocate(tmpqp)
call new(tmpqp)
call assign(tmpqp,qprop)
call fft_qprop(tmpqp%prop)
call convolution(snk_wf(wf_table(iqk,2))%wavefunc,tmpqp%prop)
call convolution(snk_wf(wf_table(iqk,2))%wavefunc,tmpqp%prop) ! double smearing
call invfft_qprop(tmpqp%prop)
call get_meson_dloop(meson_disc(isrc),Z2_noise_wf(isrc),tmpqp)
call save_meson_dloop(meson_disc(isrc),inoise,ifm,iwf)
ifm=ifm-1
call delete(tmpqp)
deallocate(tmpqp)
if (nodeid==0) write(*,'("End Disconnected Quark Trace Loop")')
if (nodeid==0) then
write(*,'("End ==== Quark #",I3," Z2 noise #",I3,A)')iqk,inoise,TRIM(str(isrc))
write(*,'("")')
write(*,*)
endif
call delete(qprop)
enddo ! end of do isrc
ifm=ifm+2
enddo ! end of do iqk
do isrc=1,2
call delete_meson_dloop(meson_disc(isrc))
call delete(Z2_noise_wf(isrc))
enddo
enddo ! end noise
if (nodeid==0) then
write(*,'("End Disconnected part")')
endif
do isrc=1,3
call delete(snk_wf(isrc))
enddo
deallocate(Z2_noise_wf)
deallocate(snk_wf)
enddo ! end of do inoise
call delete(qprop)
deallocate(qprop)
call save(main)
call delete(main)
stop
end program