2pt with source at tfl

src/io.jl

@@ -146,6 +146,8 @@ function save_data()
         for fl in 1:length(flow_times)
             BDIO_write!(fb,pp_corr_t[:,noi,fl])
             BDIO_write!(fb,ap_corr_t[:,noi,fl])
+            BDIO_write!(fb,pp_corr_tfl[:,noi,fl])
+            BDIO_write!(fb,ap_corr_tfl[:,noi,fl])
             BDIO_write!(fb,Quark_cond[:,noi,fl])
             BDIO_write!(fb,Quark_cond_cfl[:,noi,fl])
             BDIO_write!(fb,[phat_t[noi,fl]])

src/meas.jl

@@ -21,6 +21,9 @@ function load_fields()
     global pp_corr_t = fill(zero(Float64),(lp.iL[4],N_noise,length(flow_times)));
     global ap_corr_t = fill(zero(ComplexF64),(lp.iL[4],N_noise,length(flow_times)));
 
+    global pp_corr_tfl = fill(zero(Float64),(lp.iL[4],N_noise,length(flow_times)));
+    global ap_corr_tfl = fill(zero(ComplexF64),(lp.iL[4],N_noise,length(flow_times)));
+
     global Eoft = Array{Complex{Float64}}(undef,1+length(flow_times),lp.iL[4]);
     global Eofpla = Array{Complex{Float64}}(undef,lp.iL[4],lp.npls);
 
@@ -108,8 +111,8 @@ end
 """
 function one_pt()
 
-Stores the one point function in the variable 'Quark_cond' and 'Quark_cond_cfl'.
-The final result is 'Quark_cond' + C_fl*'Quark_cond_cfl'
+Stores the one point function in the variable 'Quark_cond' and 'Quark_cond_cfl'. The final result for the quark condensate is 'Quark_cond' + C_fl*'Quark_cond_cfl'.
+The two point function with source in possitive flow time is also stored in 'pp_corr_tfl' and 'ap_corr_tfl'
 """
 function one_pt()
 
@@ -149,6 +152,20 @@ function one_pt()
                     CUDA.@allowscalar Quark_cond_cfl[t,noi,fl] += -norm2(dws.sp[b,r])
                 end end end
             end
+
+            pp_density .= Array(norm2.(psi))
+            ap_density .= Array(dot.(psi,dmul.(Gamma{4},psi)))
+            pp_corr_tfl[:,noi,fl] .= zero(Float64)
+            ap_corr_tfl[:,noi,fl] .= zero(ComplexF64)
+            for t in 1:lp.iL[4]
+                for i in 1:lp.iL[1] for j in 1:lp.iL[2] for k in 1:lp.iL[3]
+                    b,r = point_index(CartesianIndex{lp.ndim}((i,j,k,t)),lp)
+                    pp_corr_tfl[t,noi,fl] += pp_density[b,r]
+                    ap_corr_tfl[t,noi,fl] += ap_density[b,r]
+                end end end
+            end
+
+
             @timeit "CPU to GPU" copyto!(U,U_CPU)
 
         end
@@ -156,6 +173,10 @@ function one_pt()
 
     Quark_cond .= Quark_cond./prod(lp.iL[1:3])
     Quark_cond_cfl .= Quark_cond_cfl./prod(lp.iL[1:3])
+
+    pp_corr_tfl .= pp_corr_tfl./prod(lp.iL[1:3])
+    ap_corr_tfl .= ap_corr_tfl./prod(lp.iL[1:3])
+
     println(log_file,"Finished measuring 1pt")
 
     return nothing

src/utils.jl

@@ -47,6 +47,8 @@ function read_ff(name::String)
     ap_corr_mc = Array{ComplexF64}(undef,0,iL[4],N_noise);
     pp_corr_t_mc = Array{Float64}(undef,0,iL[4],N_noise,nflow);
     ap_corr_t_mc = Array{ComplexF64}(undef,0,iL[4],N_noise,nflow);
+    pp_corr_tfl_mc = Array{Float64}(undef,0,iL[4],N_noise,nflow);
+    ap_corr_tfl_mc = Array{ComplexF64}(undef,0,iL[4],N_noise,nflow);
     Sigma_mc = Array{ComplexF64}(undef,0,iL[4],N_noise,nflow);
     Phat_t_mc = Array{Float64}(undef,0,N_noise,nflow);
     Sigma_cfl_mc = Array{ComplexF64}(undef,0,iL[4],N_noise,nflow);
@@ -72,6 +74,12 @@ function read_ff(name::String)
                 BDIO_read(file,TvecC)
                 ap_corr_t[1,:,noi,fl].= TvecC
 
+                BDIO_read(file,TvecR)
+                pp_corr_tfl[1,:,noi,fl].= TvecR
+
+                BDIO_read(file,TvecC)
+                ap_corr_tfl[1,:,noi,fl].= TvecC
+
                 BDIO_read(file,TvecC)
                 Sigma[1,:,noi,fl] .= TvecC
                 BDIO_read(file,TvecC)
@@ -90,6 +98,8 @@ function read_ff(name::String)
         ap_corr_mc=[ap_corr_mc;ap_corr]
         pp_corr_t_mc=[pp_corr_t_mc;pp_corr_t]
         ap_corr_t_mc=[ap_corr_t_mc;ap_corr_t]
+        pp_corr_tfl_mc=[pp_corr_tfl_mc;pp_corr_tfl]
+        ap_corr_tfl_mc=[ap_corr_tfl_mc;ap_corr_tfl]
         Sigma_mc=[Sigma_mc;Sigma]
         Sigma_cfl_mc=[Sigma_cfl_mc;Sigma_cfl]
         Phat_t_mc=[Phat_t_mc;Phat_t]
@@ -100,7 +110,7 @@ function read_ff(name::String)
 
     BDIO_close!(file)
 
-    return Eoft_mc,pp_corr_mc,ap_corr_mc,pp_corr_t_mc,ap_corr_t_mc,Sigma_mc,Sigma_cfl_mc,Phat_t_mc
+    return Eoft_mc,pp_corr_mc,ap_corr_mc,pp_corr_t_mc,ap_corr_t_mc,pp_corr_tfl_mc,ap_corr_tfl_mc,Sigma_mc,Sigma_cfl_mc,Phat_t_mc
 end
 
 
@@ -111,7 +121,7 @@ Loads in the variables 'Eoft', 'pp_corr', 'ap_corr', 'pp_corr_t','ap_corr_t', 'S
 'uwreal' quantities. The indices are, when corresponding, eucliden time and flow times
 """
 function uwff(file::String)
-    Eoft_mc,pp_corr_mc,ap_corr_mc,pp_corr_t_mc,ap_corr_t_mc,Sigma_mc,Sigma_cfl_mc,Phat_t_mc = read_ff(file)
+    Eoft_mc,pp_corr_mc,ap_corr_mc,pp_corr_t_mc,ap_corr_t_mc,pp_corr_tfl_mc,ap_corr_tfl_mc,Sigma_mc,Sigma_cfl_mc,Phat_t_mc = read_ff(file)
     runame = String(split(split(file,"/")[end],".")[1])
 
     T = size(pp_corr_mc)[2]
@@ -126,6 +136,9 @@ function uwff(file::String)
     global pp_corr_t = [uwreal(sum(pp_corr_t_mc,dims = 3)[:,i,1,k]./Ns,runame) for i in 1:T, k in 1:Nfl]
     global ap_corr_t = [uwreal(real.(sum(ap_corr_t_mc,dims = 3)[:,i,1,k]./Ns),runame) for i in 1:T,k in 1:Nfl]
 
+    global pp_corr_tfl = [uwreal(sum(pp_corr_tfl_mc,dims = 3)[:,i,1,k]./Ns,runame) for i in 1:T, k in 1:Nfl]
+    global ap_corr_tfl = [uwreal(real.(sum(ap_corr_tfl_mc,dims = 3)[:,i,1,k]./Ns),runame) for i in 1:T,k in 1:Nfl]
+
     global Sigma = [uwreal(real.(sum(Sigma_mc,dims = 3)[:,i,1,k]./Ns),runame) for i in 1:T, k in 1:Nfl]
     global Sigma_cfl = [uwreal(real.(sum(Sigma_cfl_mc,dims = 3)[:,i,1,k]./Ns),runame) for i in 1:T, k in 1:Nfl]