Fist step to split front and backflow

parent 2a36a27f
...@@ -22,11 +22,18 @@ theta = 0.0 ...@@ -22,11 +22,18 @@ theta = 0.0
theta_t = 3.1415926535897 theta_t = 3.1415926535897
csw = 1.76923076923077 csw = 1.76923076923077
ct = 9.8718528389891003e-01 ct = 9.8718528389891003e-01
tolerance = 1.0e-13
maxiter = 50000
[Frontflow]
N_noise = 6
t_zero = 2.0
epsilon = 0.01
nsteps = 100
tsource = 8 tsource = 8
[Measurements] [Backflow]
N_noise = 6 N_noise = 6
Flow_times = [2.47,3.86] Flow_times = [2.47,3.86]
Nsaves = 25 Nsaves = 25
tolerance = 1.0e-13 tsource = 8
maxiter = 50000
...@@ -63,8 +63,7 @@ function load_structs() ...@@ -63,8 +63,7 @@ function load_structs()
global int = wfl_rk3(Float64, 0.1, 1.0) global int = wfl_rk3(Float64, 0.1, 1.0)
global intsch = omf4(Float64,params["HMC"]["eps"], params["HMC"]["ns"]); global intsch = omf4(Float64,params["HMC"]["eps"], params["HMC"]["ns"]);
global flow_times = params["Measurements"]["Flow_times"] global flow_times = params["Backflow"]["Flow_times"]
global N_noise = params["Measurements"]["N_noise"]
global MCid = 0 global MCid = 0
if params["Space"]["bc"] == 0 if params["Space"]["bc"] == 0
global V4 = prod(lp.iL[1:3])*(lp.iL[4]-1) global V4 = prod(lp.iL[1:3])*(lp.iL[4]-1)
...@@ -97,11 +96,14 @@ function write_log() ...@@ -97,11 +96,14 @@ function write_log()
println(log_file,"kappa = ", params["Fermion"]["kappa"]) println(log_file,"kappa = ", params["Fermion"]["kappa"])
println(log_file,"theta = ", params["Fermion"]["theta"]) println(log_file,"theta = ", params["Fermion"]["theta"])
println(log_file,"csw = ", dpar.csw) println(log_file,"csw = ", dpar.csw)
println(log_file,"tolerance = ", params["Measurements"]["tolerance"]) println(log_file,"tolerance = ", params["Fermion"]["tolerance"])
println(log_file,"maxiter = ", params["Measurements"]["maxiter"]) println(log_file,"maxiter = ", params["Fermion"]["maxiter"])
println(log_file,"Flow times = ", flow_times) println(log_file,"N noise frontflow = ", params["Frontflow"]["N_noise"])
println(log_file,"N noise = ", N_noise) println(log_file,"eps frontflw = ", params["Frontflow"]["epsilon"])
println(log_file,"t_source = ", params["Fermion"]["tsource"]) println(log_file,"t_source frontflow = ", params["Frontflow"]["tsource"])
println(log_file,"Flow times = ", params["Backflow"]["Flow_times"])
println(log_file,"N noise backflow = ", params["Backflow"]["N_noise"])
println(log_file,"t_source backflow = ", params["Backflow"]["tsource"])
flush(log_file) flush(log_file)
return nothing return nothing
...@@ -129,9 +131,11 @@ function save_data() ...@@ -129,9 +131,11 @@ function save_data()
BDIO_write!(fb, [convert(Int32, lp.iL[i]) for i in 1:4]) BDIO_write!(fb, [convert(Int32, lp.iL[i]) for i in 1:4])
BDIO_write!(fb, [dpar.m0, dpar.csw]) BDIO_write!(fb, [dpar.m0, dpar.csw])
BDIO_write!(fb, [dpar.th[i] for i in 1:4]) BDIO_write!(fb, [dpar.th[i] for i in 1:4])
BDIO_write!(fb, [convert(Int32,params["Frontflow"]["N_noise"])])
BDIO_write!(fb, [convert(Int32,params["Frontflow"]["nsteps"])])
BDIO_write!(fb, [convert(Int32,params["Backflow"]["N_noise"])])
BDIO_write!(fb, [convert(Int32,length(flow_times))]) BDIO_write!(fb, [convert(Int32,length(flow_times))])
length(flow_times) > 0 ? BDIO_write!(fb, flow_times) : nothing length(flow_times) > 0 ? BDIO_write!(fb, flow_times) : nothing
BDIO_write!(fb, [convert(Int32,N_noise)])
BDIO_write_hash!(fb) BDIO_write_hash!(fb)
end end
...@@ -139,14 +143,23 @@ function save_data() ...@@ -139,14 +143,23 @@ function save_data()
BDIO_start_record!(fb, BDIO_BIN_GENERIC, 8) BDIO_start_record!(fb, BDIO_BIN_GENERIC, 8)
for noi in 1:N_noise for noi in 1:params["Frontflow"]["N_noise"]
BDIO_write!(fb,pp_corr_t0[:,noi]) BDIO_write!(fb,pp_corr_t0[:,noi])
BDIO_write!(fb,ap_corr_t0[:,noi]) BDIO_write!(fb,ap_corr_t0[:,noi])
for fl in 1:length(flow_times) for fl in 1:params["Frontflow"]["nsteps"]+1
BDIO_write!(fb,pp_corr_t[:,noi,fl]) BDIO_write!(fb,pp_corr_t[:,noi,fl])
BDIO_write!(fb,ap_corr_t[:,noi,fl]) BDIO_write!(fb,ap_corr_t[:,noi,fl])
BDIO_write!(fb,[phat_t[noi,fl]])
end
end
BDIO_write!(fb,Eoft)
for noi in 1:params["Backflow"]["N_noise"]
for fl in 1:length(flow_times)
BDIO_write!(fb,pp_corr_tfl[:,noi,fl]) BDIO_write!(fb,pp_corr_tfl[:,noi,fl])
BDIO_write!(fb,ap_corr_tfl[:,noi,fl]) BDIO_write!(fb,ap_corr_tfl[:,noi,fl])
...@@ -162,7 +175,6 @@ function save_data() ...@@ -162,7 +175,6 @@ function save_data()
BDIO_write!(fb,[phat_t[noi,fl]]) BDIO_write!(fb,[phat_t[noi,fl]])
end end
end end
BDIO_write!(fb,Eoft)
BDIO_write_hash!(fb) BDIO_write_hash!(fb)
......
...@@ -16,27 +16,26 @@ function load_fields() ...@@ -16,27 +16,26 @@ function load_fields()
global pp_density = Array{Float64}(undef,lp.bsz,lp.rsz); global pp_density = Array{Float64}(undef,lp.bsz,lp.rsz);
global ap_density = Array{Complex{Float64}}(undef,lp.bsz,lp.rsz); global ap_density = Array{Complex{Float64}}(undef,lp.bsz,lp.rsz);
global pp_corr_t0 = fill(zero(Float64),(lp.iL[4],N_noise)); global pp_corr_t0 = fill(zero(Float64),(lp.iL[4],params["Frontflow"]["N_noise"]));
global ap_corr_t0 = fill(zero(ComplexF64),(lp.iL[4],N_noise)); global ap_corr_t0 = fill(zero(ComplexF64),(lp.iL[4],params["Frontflow"]["N_noise"]));
global pp_corr_t = fill(zero(Float64),(lp.iL[4],N_noise,length(flow_times))); global pp_corr_t = fill(zero(Float64),(lp.iL[4],params["Frontflow"]["N_noise"],params["Frontflow"]["nsteps"]+1));
global ap_corr_t = fill(zero(ComplexF64),(lp.iL[4],N_noise,length(flow_times))); global ap_corr_t = fill(zero(ComplexF64),(lp.iL[4],params["Frontflow"]["N_noise"],params["Frontflow"]["nsteps"]+1));
global pp_corr_tfl = fill(zero(Float64),(lp.iL[4],N_noise,length(flow_times))); global pp_corr_tfl = fill(zero(Float64),(lp.iL[4],params["Backflow"]["N_noise"],length(flow_times)));
global ap_corr_tfl = fill(zero(ComplexF64),(lp.iL[4],N_noise,length(flow_times))); global ap_corr_tfl = fill(zero(ComplexF64),(lp.iL[4],params["Backflow"]["N_noise"],length(flow_times)));
global Eoft = Array{Complex{Float64}}(undef,1+length(flow_times),lp.iL[4]); global Eoft = Array{Complex{Float64}}(undef,2+ffnsteps,lp.iL[4]);
global Eofpla = Array{Complex{Float64}}(undef,lp.iL[4],lp.npls); global Eofpla = Array{Complex{Float64}}(undef,lp.iL[4],lp.npls);
global phat_t = Array{Float64}(undef,N_noise,length(flow_times)); global phat_t = Array{Float64}(undef,params["Frontflow"]["N_noise"],length(flow_times));
global Quark_cond = Array{Complex{Float64}}(undef,lp.iL[4],N_noise,length(flow_times)); global Quark_cond = Array{Complex{Float64}}(undef,lp.iL[4],params["Backflow"]["N_noise"],length(flow_times));
global Quark_cond_cfl = Array{Complex{Float64}}(undef,lp.iL[4],N_noise,length(flow_times)); global Quark_cond_cfl = Array{Complex{Float64}}(undef,lp.iL[4],params["Backflow"]["N_noise"],length(flow_times));
global Quark_cond2 = Array{Complex{Float64}}(undef,lp.iL[4],params["Backflow"]["N_noise"],length(flow_times));
global Quark_cond2 = Array{Complex{Float64}}(undef,lp.iL[4],N_noise,length(flow_times)); global Quark_cond2_cfl = Array{Complex{Float64}}(undef,lp.iL[4],params["Backflow"]["N_noise"],length(flow_times));
global Quark_cond2_cfl = Array{Complex{Float64}}(undef,lp.iL[4],N_noise,length(flow_times)); global ChiDchi = Array{Complex{Float64}}(undef,lp.iL[4],params["Backflow"]["N_noise"],length(flow_times));
global ChiDchi = Array{Complex{Float64}}(undef,lp.iL[4],N_noise,length(flow_times)); global ChiDchi_cfl = Array{Complex{Float64}}(undef,lp.iL[4],params["Backflow"]["N_noise"],length(flow_times));
global ChiDchi_cfl = Array{Complex{Float64}}(undef,lp.iL[4],N_noise,length(flow_times));
return nothing return nothing
end end
...@@ -46,17 +45,17 @@ function two_pt() ...@@ -46,17 +45,17 @@ function two_pt()
Stores the two point function in the variables 'pp_corr_t0', 'ap_corr_t0', 'pp_corr_t' and 'ap_corr_t' Stores the two point function in the variables 'pp_corr_t0', 'ap_corr_t0', 'pp_corr_t' and 'ap_corr_t'
""" """
function two_pt() function two_pt() # Will be Frontflow
println(log_file,"Measuring 2pt...") println(log_file,"\nMeasuring 2pt...")
flush(log_file) flush(log_file)
Eoft_plaq(Eofpla,U, gp, lp, ymws) Eoft_plaq(Eofpla,U, gp, lp, ymws)
Eoft[1,:] .= sum(Eofpla,dims = 2) Eoft[1,:] .= sum(Eofpla,dims = 2)
for noi in 1:N_noise for noi in 1:params["Frontflow"]["N_noise"]
niter = propagator!(psi,U, dpar, dws, lp,params["Measurements"]["maxiter"], params["Measurements"]["tolerance"],params["Fermion"]["tsource"]) niter = propagator!(psi,U, dpar, dws, lp,params["Fermion"]["maxiter"], params["Fermion"]["tolerance"],params["Frontflow"]["tsource"])
println(log_file,"CG converged after ",niter," iterations with residue ",CUDA.mapreduce(x -> norm2(x), +, dws.sr)) println(log_file,"CG converged after ",niter," iterations with residue ",CUDA.mapreduce(x -> norm2(x), +, dws.sr))
flush(log_file) flush(log_file)
...@@ -74,28 +73,40 @@ function two_pt() ...@@ -74,28 +73,40 @@ function two_pt()
end end end end end end
end end
delta_t = flow_times .- [0.0,flow_times[1:end-1]...] _,epslist = flw_adapt(U, psi, int, params["Frontflow"]["t_zero"], gp, dpar, lp, ymws, dws)
neps = int.eps_ini
for fl in 1:length(flow_times) for tstep in 1:params["Frontflow"]["nsteps"]
_,epslist = flw_adapt(U, psi, int, delta_t[fl], neps, gp, dpar, lp, ymws, dws)
neps = epslist[end]
if noi ==1 if noi ==1
Eoft_plaq(Eofpla,U, gp, lp, ymws) Eoft_plaq(Eofpla,U, gp, lp, ymws)
Eoft[fl+1,:] .= sum(Eofpla,dims = 2) Eoft[tstep+1,:] .= sum(Eofpla,dims = 2)
end end
pp_density .= Array(norm2.(psi)) pp_density .= Array(norm2.(psi))
ap_density .= Array(dot.(psi,dmul.(Gamma{4},psi))) ap_density .= Array(dot.(psi,dmul.(Gamma{4},psi)))
pp_corr_t[:,noi,fl] .= zero(Float64) pp_corr_t[:,noi,tstep] .= zero(Float64)
ap_corr_t[:,noi,fl] .= zero(ComplexF64) ap_corr_t[:,noi,tstep] .= zero(ComplexF64)
for t in 1:lp.iL[4] 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] 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) b,r = point_index(CartesianIndex{lp.ndim}((i,j,k,t)),lp)
pp_corr_t[t,noi,fl] += pp_density[b,r] pp_corr_t[t,noi,tstep] += pp_density[b,r]
ap_corr_t[t,noi,fl] += ap_density[b,r] ap_corr_t[t,noi,tstep] += ap_density[b,r]
end end end end end end
end end
flw(U, psi, int, 1, params["Frontflow"]["epsilon"], gp, dpar, lp, ymws, dws)
end
pp_density .= Array(norm2.(psi))
ap_density .= Array(dot.(psi,dmul.(Gamma{4},psi)))
pp_corr_t[:,noi,end] .= zero(Float64)
ap_corr_t[:,noi,end] .= 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_t[t,noi,end] += pp_density[b,r]
ap_corr_t[t,noi,end] += ap_density[b,r]
end end end
end end
@timeit "CPU to GPU" copyto!(U,U_CPU) @timeit "CPU to GPU" copyto!(U,U_CPU)
...@@ -118,22 +129,22 @@ function one_pt() ...@@ -118,22 +129,22 @@ function one_pt()
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'. 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' The two point function with source in possitive flow time is also stored in 'pp_corr_tfl' and 'ap_corr_tfl'
""" """
function one_pt() function one_pt() # Will be backflow
@timeit "CPU to GPU" copyto!(U,U_CPU) @timeit "CPU to GPU" copyto!(U,U_CPU)
for fl in 1:length(flow_times) for fl in 1:length(flow_times)
println(log_file,"Measuring 1pt at t = "*string(flow_times[fl])*"...") println(log_file,"\nMeasuring 1pt at t = "*string(flow_times[fl])*"...")
flush(log_file) flush(log_file)
for noi in 1:N_noise for noi in 1:params["Backflow"]["N_noise"]
## The change in the field before randomize should be moved to the package ## The change in the field before randomize should be moved to the package
psi .= 0.0*psi psi .= 0.0*psi
pfrandomize!(psi,lp,params["Fermion"]["tsource"]) pfrandomize!(psi,lp,params["Backflow"]["tsource"])
@timeit "GPU to CPU" psit_CPU .= Array(psi) @timeit "GPU to CPU" psit_CPU .= Array(psi)
backflow(psi, U, flow_times[fl], params["Measurements"]["Nsaves"], gp, dpar, lp, ymws, dws) backflow(psi, U, flow_times[fl], , gp, dpar, lp, ymws, dws)
@timeit "GPU to CPU" psi_CPU .= Array(psi) @timeit "GPU to CPU" psi_CPU .= Array(psi)
...@@ -141,7 +152,7 @@ function one_pt() ...@@ -141,7 +152,7 @@ function one_pt()
dws.sp .= dmul.(Gamma{5},psi) dws.sp .= dmul.(Gamma{5},psi)
g5Dw!(psi,U,dws.sp,dpar,dws,lp) g5Dw!(psi,U,dws.sp,dpar,dws,lp)
niter = CG!(psi,U,DwdagDw!,dpar,lp,dws,params["Measurements"]["maxiter"], params["Measurements"]["tolerance"]) niter = CG!(psi,U,DwdagDw!,dpar,lp,dws,params["Fermion"]["maxiter"], params["Fermion"]["tolerance"])
end end
println(log_file,"CG converged after ",niter," iterations with residue ",CUDA.mapreduce(x -> norm2(x), +, dws.sr)) println(log_file,"CG converged after ",niter," iterations with residue ",CUDA.mapreduce(x -> norm2(x), +, dws.sr))
...@@ -205,9 +216,6 @@ function one_pt() ...@@ -205,9 +216,6 @@ function one_pt()
end end end end end end
end end
Quark_cond2[:,noi,fl] .= zero(ComplexF64) Quark_cond2[:,noi,fl] .= zero(ComplexF64)
ap_density .= Array(dot(dws.st[b,r],psi[b,r])) ap_density .= Array(dot(dws.st[b,r],psi[b,r]))
for t in 1:lp.iL[4] for t in 1:lp.iL[4]
...@@ -285,7 +293,7 @@ function two_pt_rnd() ...@@ -285,7 +293,7 @@ function two_pt_rnd()
println(log_file,"Measuring PhatPt correlator...") println(log_file,"Measuring PhatPt correlator...")
flush(log_file) flush(log_file)
for noi in 1:N_noise for noi in 1:params["Frontflow"]["N_noise"]
pfrandomize!(psi,lp) pfrandomize!(psi,lp)
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment