Main script and dirs. Only four correlation functions.

This commit is contained in:
Fernando Pérez Panadero 2024-02-19 17:20:28 +01:00
commit 713812c3b9
2 changed files with 232 additions and 0 deletions

21
input/sfcf.in Normal file
View file

@ -0,0 +1,21 @@
[Run]
user = "fperez"
name = "SFtest_Setup10"
[Space]
size = [8,8,8,8]
blocks = [4,4,4,4]
phi0 = [-1.047197551196598, 0.0] # Boundary fields for SF
phiT = [-3.141592653589793, 1.047197551196598]
cG = 1.0 # Correction for boundaries (SF only)
[Fermion]
beta = 6.0 #Used only if cG/csw/ct are not input
kappa = 0.1
theta = 0.5
csw = 1.0
ct = 0.9
[Solver]
tolerance = 1.0e-10
maxiter = 1000

211
sfcf.jl Normal file
View file

@ -0,0 +1,211 @@
using LatticeGPU
using TOML
using TimerOutputs
using ArgParse
using CUDA
function parse_commandline()
s = ArgParseSettings()
@add_arg_table s begin
"-i"
help = "Input parameters file"
required = true
arg_type = String
"-c"
help = "Gauge configuration file"
required = true
arg_type = String
end
return parse_args(s)
end
parsed_args = parse_commandline()
println("--------------------------------------------------")
println("Reading input file from:", parsed_args["i"], "...")
params = TOML.parsefile(parsed_args["i"])
cG = try
params["Space"]["cG"]
catch
try
println("Using default value for cG")
1.0 - 0.089*6/params["Fermion"]["beta"]
catch
error("Either cG of beta has to be defined in the input file")
end
end
ct = try
params["Fermion"]["ct"]
catch
try
error("ct not defined")
println("Using default value for ct")
catch
error("Either ct of beta has to be defined in the input file")
end
end
csw = try
params["Fermion"]["csw"]
catch
try
error("csw not defined")
println("Using default value for csw")
catch
error("Either csw of beta has to be defined in the input file")
end
end
lp = SpaceParm{4}(tuple(params["Space"]["size"]...), tuple(params["Space"]["blocks"]...),BC_SF_ORBI, (0,0,0,0,0,0))
gp = GaugeParm{Float64}(SU3{Float64},params["Fermion"]["beta"],1.0,(cG,0.0),params["Space"]["phiT"],lp.iL);
dpar = DiracParam{Float64}(SU3fund,(1/(2*params["Fermion"]["kappa"])) - 4,csw,ntuple(i -> exp((i!=4)*im*params["Fermion"]["theta"]/lp.iL[i]),4),0.0,ct);
dws = DiracWorkspace(SU3fund,Float64,lp);
ymws = YMworkspace(SU3,Float64,lp);
println("Reading gauge field from: ", parsed_args["c"], "...")
U,_ = read_cnfg(parsed_args["c"])
Csw!(dws, U, gp, lp)
println("--------------------------------------------------")
println("Parameters:")
println("Lattice size: ", lp.iL)
println("Phi0 = ", params["Space"]["phi0"])
println("PhiT = ", params["Space"]["phiT"])
println("cG = ", gp.cG[1])
println("kappa = ", params["Fermion"]["kappa"])
println("theta = ", params["Fermion"]["theta"])
println("csw = ", dpar.csw)
println("ct = ", dpar.ct)
println("tolerance = ", params["Solver"]["tolerance"])
println("maxiter = ", params["Solver"]["maxiter"])
println("--------------------------------------------------")
psi = scalar_field(Spinor{4,SU3fund{Float64}},lp)
f1 = 0.0
println("Computing propagators from t=0 to the bulk")
niter = bndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 1, 1)
println("CG converged for s=1 and c=1 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
A11 = Array(psi)
f1 += norm2(bndtobnd(psi, U, dpar, dws, lp))
niter = bndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 2, 1)
println("CG converged for c=2 and s=1 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
A21 = Array(psi)
f1 += norm2(bndtobnd(psi, U, dpar, dws, lp))
niter = bndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 3, 1)
println("CG converged for c=3 and s=1 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
A31 = Array(psi)
f1 += norm2(bndtobnd(psi, U, dpar, dws, lp))
niter = bndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 1, 2)
println("CG converged for c=1 and s=2 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
A12 = Array(psi)
f1 += norm2(bndtobnd(psi, U, dpar, dws, lp))
niter = bndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 2, 2)
println("CG converged for c=2 and s=2 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
A22 = Array(psi)
f1 += norm2(bndtobnd(psi, U, dpar, dws, lp))
niter = bndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 3, 2)
println("CG converged for c=3 and s=2 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
A32 = Array(psi)
f1 += norm2(bndtobnd(psi, U, dpar, dws, lp))
println("Computing propagators from t=T to the bulk")
niter = Tbndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 1, 1)
println("CG converged for s=1 and c=1 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
B11 = Array(psi)
niter = Tbndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 2, 1)
println("CG converged for c=2 and s=1 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
B21 = Array(psi)
niter = Tbndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 3, 1)
println("CG converged for c=3 and s=1 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
B31 = Array(psi)
niter = Tbndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 1, 2)
println("CG converged for c=1 and s=2 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
B12 = Array(psi)
niter = Tbndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 2, 2)
println("CG converged for c=2 and s=2 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
B22 = Array(psi)
niter = Tbndpropagator!(psi, U, dpar, dws, lp, params["Solver"]["maxiter"], params["Solver"]["tolerance"], 3, 2)
println("CG converged for c=3 and s=2 after ",niter," iterations with absolute residue ", CUDA.mapreduce(x -> norm(x), +, dws.sr))
B32 = Array(psi)
println("--------------------------------------------------")
println("Computing correlators")
println("Computing fP...")
file = open(basename(parsed_args["i"])*".fP","w+")
fP = Vector{Complex{Float64}}(undef,lp.iL[4])
for t in 1:lp.iL[4]
fP[t] = 0.0
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)
fP[t] += norm2(A11[b,r]) + norm2(A21[b,r]) + norm2(A31[b,r]) + norm2(A12[b,r]) + norm2(A22[b,r]) + norm2(A32[b,r])
end end end
println(file,t," ",fP[t])
end
flush(file)
close(file)
println("Computing fA...")
file = open(basename(parsed_args["i"])*".fA","w+")
fA = Vector{Complex{Float64}}(undef,lp.iL[4])
for t in 1:lp.iL[4]
fA[t] = 0.0
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)
fA[t] += dot(A11[b,r],A13[b,r]) + dot(A21[b,r],A23[b,r]) + dot(A31[b,r],A33[b,r]) + dot(A12[b,r],A14[b,r]) + dot(A22[b,r],A24[b,r]) + dot(A32[b,r],A34[b,r])
end end end
println(file,t," ",fA[t])
end
flush(file)
close(file)
println("Computing gP...")
file = open(basename(parsed_args["i"])*".gP","w+")
gP = Vector{Complex{Float64}}(undef,lp.iL[4])
for t in 1:lp.iL[4]
gP[t] = 0.0
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)
gP[t] += norm2(B11[b,r]) + norm2(B21[b,r]) + norm2(B31[b,r]) + norm2(B12[b,r]) + norm2(B22[b,r]) + norm2(B32[b,r])
println(file,t," ",gP[t])
end end end
end
flush(file)
close(file)
println("Computing gA...")
file = open(basename(parsed_args["i"])*".gA","w+")
gA = Vector{Complex{Float64}}(undef,lp.iL[4])
for t in 1:lp.iL[4]
gA[t] = 0.0
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)
gA[t] += dot(B11[b,r],B13[b,r]) + dot(B21[b,r],B23[b,r]) + dot(B31[b,r],B33[b,r]) + dot(B12[b,r],B14[b,r]) + dot(B22[b,r],B24[b,r]) + dot(B32[b,r],B34[b,r])
end end end
println(file,t," ",gA[t])
end
flush(file)
close(file)