mirror of
https://igit.ific.uv.es/alramos/latticegpu.jl.git
synced 2025-09-18 12:17:14 +02:00
Even/Odd splitting in OR sweeps
This commit is contained in:
Alberto Ramos
parent
eabbfe6307
commit
bdb69f5d0e
4 changed files with 222 additions and 48 deletions
|
|
@ -47,7 +47,7 @@ export FlowIntr, wfl_euler, zfl_euler, wfl_rk2, zfl_rk2, wfl_rk3, zfl_rk3
|
||||||
export flw, flw_adapt
|
export flw, flw_adapt
|
||||||
export sfcoupling, bndfield, setbndfield
|
export sfcoupling, bndfield, setbndfield
|
||||||
export import_lex64, import_cern64, import_bsfqcd, save_cnfg, read_cnfg, read_gp
|
export import_lex64, import_cern64, import_bsfqcd, save_cnfg, read_cnfg, read_gp
|
||||||
export updt_or_wilson!
|
export evenodd, updt_or_wilson!
|
||||||
|
|
||||||
include("Spinors/Spinors.jl")
|
include("Spinors/Spinors.jl")
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -182,6 +182,6 @@ include("YMio.jl")
|
||||||
export import_lex64, import_cern64, import_bsfqcd, save_cnfg, read_cnfg, read_gp
|
export import_lex64, import_cern64, import_bsfqcd, save_cnfg, read_cnfg, read_gp
|
||||||
|
|
||||||
include("YMupdate.jl")
|
include("YMupdate.jl")
|
||||||
export updt_or_wilson!
|
export evenodd, updt_or_wilson!
|
||||||
|
|
||||||
end
|
end
|
||||||
|
|
|
||||||
|
|
@ -9,19 +9,83 @@
|
||||||
### created: Mon Sep 1 08:59:22 2025
|
### created: Mon Sep 1 08:59:22 2025
|
||||||
###
|
###
|
||||||
|
|
||||||
|
struct EvenOdd{T}
|
||||||
|
state::Bool
|
||||||
|
e::T
|
||||||
|
o::T
|
||||||
|
end
|
||||||
|
|
||||||
|
function evenodd(eo::EvenOdd, lp::SpaceParm)
|
||||||
|
|
||||||
|
if eo.state
|
||||||
|
return eo
|
||||||
|
end
|
||||||
|
|
||||||
|
return evenodd(lp)
|
||||||
|
end
|
||||||
|
|
||||||
|
function evenodd(lp::SpaceParm)
|
||||||
|
|
||||||
|
if all(lp.iL .% 2 .== 0)
|
||||||
|
e = Array{Int64,1}(undef, div(lp.bsz,2))
|
||||||
|
o = Array{Int64,1}(undef, div(lp.bsz,2))
|
||||||
|
|
||||||
|
ie = 1
|
||||||
|
io = 1
|
||||||
|
for b in 1:lp.bsz
|
||||||
|
clr = point_color((b, 1), lp)
|
||||||
|
if clr % 2 == 0
|
||||||
|
e[ie] = b
|
||||||
|
ie = ie + 1
|
||||||
|
else
|
||||||
|
o[io] = b
|
||||||
|
io = io + 1
|
||||||
|
end
|
||||||
|
end
|
||||||
|
return EvenOdd(true, CuArray(e), CuArray(o))
|
||||||
|
else
|
||||||
|
return EvenOdd(false,nothing, nothing)
|
||||||
|
end
|
||||||
|
|
||||||
|
end
|
||||||
|
|
||||||
|
function updt_or_wilson!(U, gp::GaugeParm, eo::EvenOdd, lp::SpaceParm{N,M,B,D}) where {N,M,B,D}
|
||||||
|
|
||||||
|
if eo.state
|
||||||
|
@timeit "OR update (Wilson action [eo])" begin
|
||||||
|
ztw = ztwist(gp, lp)
|
||||||
|
for id in 1:N
|
||||||
|
CUDA.@sync begin
|
||||||
|
CUDA.@cuda threads=div(lp.bsz,2) blocks=lp.rsz krnl_or_wilson!(U, id, eo.e, ztw, gp, lp)
|
||||||
|
end
|
||||||
|
end
|
||||||
|
|
||||||
|
for id in 1:N
|
||||||
|
CUDA.@sync begin
|
||||||
|
CUDA.@cuda threads=div(lp.bsz,2) blocks=lp.rsz krnl_or_wilson!(U, id, eo.o, ztw, gp, lp)
|
||||||
|
end
|
||||||
|
end
|
||||||
|
end
|
||||||
|
else
|
||||||
|
updt_or_wilson!(U, gp, lp)
|
||||||
|
end
|
||||||
|
|
||||||
|
return nothing
|
||||||
|
end
|
||||||
|
|
||||||
function updt_or_wilson!(U, gp::GaugeParm, lp::SpaceParm{N,M,B,D}) where {N,M,B,D}
|
function updt_or_wilson!(U, gp::GaugeParm, lp::SpaceParm{N,M,B,D}) where {N,M,B,D}
|
||||||
|
|
||||||
@timeit "OR update (Wilson action)" begin
|
@timeit "OR update (Wilson action)" begin
|
||||||
ztw = ztwist(gp, lp)
|
ztw = ztwist(gp, lp)
|
||||||
for id in 1:N
|
for id in 1:N
|
||||||
CUDA.@sync begin
|
CUDA.@sync begin
|
||||||
CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_or_wilson!(U, id, 0, ztw, lp)
|
CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_or_wilson!(U, id, 0, ztw, gp, lp)
|
||||||
end
|
end
|
||||||
end
|
end
|
||||||
|
|
||||||
for id in 1:N
|
for id in 1:N
|
||||||
CUDA.@sync begin
|
CUDA.@sync begin
|
||||||
CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_or_wilson!(U, id, 1, ztw, lp)
|
CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_or_wilson!(U, id, 1, ztw, gp, lp)
|
||||||
end
|
end
|
||||||
end
|
end
|
||||||
end
|
end
|
||||||
|
|
@ -30,19 +94,28 @@ function updt_or_wilson!(U, gp::GaugeParm, lp::SpaceParm{N,M,B,D}) where {N,M,B,
|
||||||
end
|
end
|
||||||
|
|
||||||
|
|
||||||
function krnl_or_wilson!(U::AbstractArray{SU3{T}}, id, icl, ztw, lp::SpaceParm) where {T}
|
function krnl_or_wilson!(U::AbstractArray{SU3{T}}, id, icl::Int64, ztw, gp::GaugeParm, lp::SpaceParm) where {T}
|
||||||
|
|
||||||
b = Int64(CUDA.threadIdx().x)
|
b = Int64(CUDA.threadIdx().x)
|
||||||
r = Int64(CUDA.blockIdx().x)
|
r = Int64(CUDA.blockIdx().x)
|
||||||
clr = point_color((b, r), lp)
|
clr = point_color((b, r), lp)
|
||||||
if clr % 2 == icl
|
if clr % 2 == icl
|
||||||
stp = sum_staples(U, id, (b, r), ztw, lp)
|
stp = sum_staples(U, id, (b, r), ztw, gp, lp)
|
||||||
U[b,id,r] = do_or(U[b,id,r], stp)
|
U[b,id,r] = do_or(U[b,id,r], stp)
|
||||||
end
|
end
|
||||||
|
|
||||||
return nothing
|
return nothing
|
||||||
end
|
end
|
||||||
|
|
||||||
|
function krnl_or_wilson!(U::AbstractArray{SU3{T}}, id, eo::AbstractArray, ztw, gp::GaugeParm, lp::SpaceParm) where {T}
|
||||||
|
|
||||||
|
b = eo[Int64(CUDA.threadIdx().x)]
|
||||||
|
r = Int64(CUDA.blockIdx().x)
|
||||||
|
stp = sum_staples(U, id, (b, r), ztw, gp, lp)
|
||||||
|
U[b,id,r] = do_or(U[b,id,r], stp)
|
||||||
|
|
||||||
|
return nothing
|
||||||
|
end
|
||||||
|
|
||||||
@inline function overrelax_cd(w::Tuple{T,T}) where {T}
|
@inline function overrelax_cd(w::Tuple{T,T}) where {T}
|
||||||
|
|
||||||
|
|
@ -56,44 +129,6 @@ end
|
||||||
end
|
end
|
||||||
|
|
||||||
|
|
||||||
@inline function sum_staples(U::AbstractArray{SU3{T}},im,pt::NTuple{2,Int64}, ztw, lp::SpaceParm{N,M,BC_PERIODIC,D}) where {T,N,M,D}
|
|
||||||
|
|
||||||
I = point_coord(pt, lp)
|
|
||||||
bum, rum = up(pt, im, lp)
|
|
||||||
stp = zero(M3x3{T})
|
|
||||||
for in in 1:N
|
|
||||||
if (in != im)
|
|
||||||
zf = one(ztw[begin])
|
|
||||||
for i in 1:length(lp.plidx)
|
|
||||||
if (lp.plidx[i] == (in, im)) || (lp.plidx[i] == (im, in))
|
|
||||||
if im > in
|
|
||||||
zf = ztw[i]
|
|
||||||
else
|
|
||||||
zf = conj(ztw[i])
|
|
||||||
end
|
|
||||||
break
|
|
||||||
end
|
|
||||||
end
|
|
||||||
|
|
||||||
bun, run, bdn, rdn = updw(pt, in, lp)
|
|
||||||
bd, rd = up((bdn, rdn), im, lp)
|
|
||||||
m1 = convert(M3x3{T}, U[bum,in,rum]/U[bun,im,run]/U[pt[1],in,pt[2]])
|
|
||||||
if (I[im]==1) && (I[in]==1)
|
|
||||||
m1 = m1 * zf
|
|
||||||
end
|
|
||||||
m2 = convert(M3x3{T},dag(U[bdn,in,rdn]\U[bdn,im,rdn]*U[bd, in, rd]))
|
|
||||||
Id = point_coord((bdn, rdn), lp)
|
|
||||||
if (Id[im]==1) && (Id[in]==1)
|
|
||||||
m2 = m2 * conj(zf)
|
|
||||||
end
|
|
||||||
|
|
||||||
stp = stp + m1 + m2
|
|
||||||
end
|
|
||||||
end
|
|
||||||
|
|
||||||
return stp
|
|
||||||
end
|
|
||||||
|
|
||||||
@inline function do_or(U::SU3{T}, stp) where {T}
|
@inline function do_or(U::SU3{T}, stp) where {T}
|
||||||
|
|
||||||
M = do_or12(U, stp)
|
M = do_or12(U, stp)
|
||||||
|
|
@ -172,4 +207,118 @@ end
|
||||||
return SU3{T}(u11, u12, u13, u21, u22, u23)
|
return SU3{T}(u11, u12, u13, u21, u22, u23)
|
||||||
end
|
end
|
||||||
|
|
||||||
|
@inline function sum_staples(U::AbstractArray{SU3{T}},im,pt::NTuple{2,Int64}, ztw, gp::GaugeParm, lp::SpaceParm{N,M,BC_PERIODIC,D}) where {T,N,M,D}
|
||||||
|
|
||||||
|
I = point_coord(pt, lp)
|
||||||
|
bum, rum = up(pt, im, lp)
|
||||||
|
stp = zero(M3x3{T})
|
||||||
|
for in in 1:N
|
||||||
|
if (in != im)
|
||||||
|
zf = one(ztw[begin])
|
||||||
|
for i in 1:length(lp.plidx)
|
||||||
|
if (lp.plidx[i] == (in, im)) || (lp.plidx[i] == (im, in))
|
||||||
|
if im > in
|
||||||
|
zf = ztw[i]
|
||||||
|
else
|
||||||
|
zf = conj(ztw[i])
|
||||||
|
end
|
||||||
|
break
|
||||||
|
end
|
||||||
|
end
|
||||||
|
|
||||||
|
bun, run, bdn, rdn = updw(pt, in, lp)
|
||||||
|
bd, rd = up((bdn, rdn), im, lp)
|
||||||
|
m1 = convert(M3x3{T}, U[bum,in,rum]/U[bun,im,run]/U[pt[1],in,pt[2]])
|
||||||
|
if (I[im]==1) && (I[in]==1)
|
||||||
|
m1 = m1 * zf
|
||||||
|
end
|
||||||
|
m2 = convert(M3x3{T},dag(U[bdn,in,rdn]\U[bdn,im,rdn]*U[bd, in, rd]))
|
||||||
|
Id = point_coord((bdn, rdn), lp)
|
||||||
|
if (Id[im]==1) && (Id[in]==1)
|
||||||
|
m2 = m2 * conj(zf)
|
||||||
|
end
|
||||||
|
|
||||||
|
stp = stp + m1 + m2
|
||||||
|
end
|
||||||
|
end
|
||||||
|
|
||||||
|
return stp
|
||||||
|
end
|
||||||
|
|
||||||
|
@inline function sum_staples(U::AbstractArray{SU3{T}},im,pt::NTuple{2,Int64}, ztw, gp::GaugeParm, lp::SpaceParm{N,M,BC_SF_ORBI,D}) where {T,N,M,D}
|
||||||
|
|
||||||
|
I = point_coord(pt, lp)
|
||||||
|
bum, rum = up(pt, im, lp)
|
||||||
|
stp = zero(M3x3{T})
|
||||||
|
for in in 1:N
|
||||||
|
if (in != im)
|
||||||
|
bun, run, bdn, rdn = updw(pt, in, lp)
|
||||||
|
bd, rd = up((bdn, rdn), im, lp)
|
||||||
|
|
||||||
|
if (I[end] == lp.iL[end])
|
||||||
|
if (in==N)
|
||||||
|
m1 = gp.cG[1]*convert(M3x3{T}, U[bum,in,rum]/gp.Ubnd[im]/U[pt[1],in,pt[2]])
|
||||||
|
elseif (im==N)
|
||||||
|
m1 = gp.cG[1]*convert(M3x3{T}, gp.Ubnd[in]/U[bun,im,run]/U[pt[1],in,pt[2]])
|
||||||
|
else
|
||||||
|
m1 = zero(M3x3{T})
|
||||||
|
end
|
||||||
|
else
|
||||||
|
m1 = convert(M3x3{T}, U[bum,in,rum]/U[bun,im,run]/U[pt[1],in,pt[2]])
|
||||||
|
end
|
||||||
|
|
||||||
|
if I[end] == 1
|
||||||
|
if (in==N) || (im==N)
|
||||||
|
m2 = gp.cG[2]*convert(M3x3{T},dag(U[bdn,in,rdn]\U[bdn,im,rdn]*U[bd, in, rd]))
|
||||||
|
else
|
||||||
|
m2 = zero(M3x3{T})
|
||||||
|
end
|
||||||
|
else
|
||||||
|
m2 = convert(M3x3{T},dag(U[bdn,in,rdn]\U[bdn,im,rdn]*U[bd, in, rd]))
|
||||||
|
end
|
||||||
|
stp = stp + m1 + m2
|
||||||
|
|
||||||
|
end
|
||||||
|
end
|
||||||
|
|
||||||
|
return stp
|
||||||
|
end
|
||||||
|
|
||||||
|
@inline function sum_staples(U::AbstractArray{SU3{T}},im,pt::NTuple{2,Int64}, ztw, gp::GaugeParm, lp::SpaceParm{N,M,BC_OPEN,D}) where {T,N,M,D}
|
||||||
|
|
||||||
|
I = point_coord(pt, lp)
|
||||||
|
bum, rum = up(pt, im, lp)
|
||||||
|
stp = zero(M3x3{T})
|
||||||
|
for in in 1:N
|
||||||
|
if (in != im)
|
||||||
|
bun, run, bdn, rdn = updw(pt, in, lp)
|
||||||
|
bd, rd = up((bdn, rdn), im, lp)
|
||||||
|
|
||||||
|
if (I[end] == lp.iL[end])
|
||||||
|
if (in==N)
|
||||||
|
m1 = gp.cG[1]*convert(M3x3{T}, U[bum,in,rum]/gp.Ubnd[im]/U[pt[1],in,pt[2]])
|
||||||
|
else
|
||||||
|
m1 = zero(M3x3{T})
|
||||||
|
end
|
||||||
|
else
|
||||||
|
m1 = convert(M3x3{T}, U[bum,in,rum]/U[bun,im,run]/U[pt[1],in,pt[2]])
|
||||||
|
end
|
||||||
|
|
||||||
|
if I[end] == 1
|
||||||
|
if (in==N) || (im==N)
|
||||||
|
m2 = gp.cG[2]*convert(M3x3{T},dag(U[bdn,in,rdn]\U[bdn,im,rdn]*U[bd, in, rd]))
|
||||||
|
else
|
||||||
|
m2 = zero(M3x3{T})
|
||||||
|
end
|
||||||
|
else
|
||||||
|
m2 = convert(M3x3{T},dag(U[bdn,in,rdn]\U[bdn,im,rdn]*U[bd, in, rd]))
|
||||||
|
end
|
||||||
|
stp = stp + m1 + m2
|
||||||
|
|
||||||
|
end
|
||||||
|
end
|
||||||
|
|
||||||
|
return stp
|
||||||
|
end
|
||||||
|
|
||||||
|
|
||||||
|
|
|
||||||
|
|
@ -9,26 +9,49 @@
|
||||||
### created: Tue Sep 2 16:29:48 2025
|
### created: Tue Sep 2 16:29:48 2025
|
||||||
###
|
###
|
||||||
|
|
||||||
using LatticeGPU, Test, CUDA
|
using LatticeGPU, Test, CUDA, TimerOutputs
|
||||||
|
|
||||||
T = Float64
|
T = Float64
|
||||||
lp = SpaceParm{4}((4,4,4,4), (4,4,4,4), BC_PERIODIC, (0,0,0,0,0,1))
|
#lp = SpaceParm{4}((4,4,4,4), (4,4,4,4), BC_SF_ORBI, (0,0,0,0,0,0))
|
||||||
gp = GaugeParm{T}(SU3{T}, 6.5, 1.0)
|
lp = SpaceParm{4}((16,16,16,16), (4,4,4,4), BC_PERIODIC, (0,0,0,0,0,1))
|
||||||
|
gp = GaugeParm{T}(SU3{T}, 6.5, 1.0, (1.5,1.5), (0.0,0.0), lp.iL[1])
|
||||||
ymws = YMworkspace(SU3, T, lp)
|
ymws = YMworkspace(SU3, T, lp)
|
||||||
|
|
||||||
|
println(gp.Ubnd)
|
||||||
randomize!(ymws.mom, lp, ymws)
|
randomize!(ymws.mom, lp, ymws)
|
||||||
U = exp.(ymws.mom)
|
U = exp.(ymws.mom)
|
||||||
act = gauge_action(U, lp, gp, ymws)
|
act = gauge_action(U, lp, gp, ymws)
|
||||||
plq = plaquette(U, lp, gp, ymws)
|
plq = plaquette(U, lp, gp, ymws)
|
||||||
pl_exact = Eoft_plaq(U, gp, lp, ymws)
|
pl_exact = Eoft_plaq(U, gp, lp, ymws)
|
||||||
cl_exact = Eoft_clover(U, gp, lp, ymws)
|
cl_exact = Eoft_clover(U, gp, lp, ymws)
|
||||||
|
println(lp)
|
||||||
println("## Random config: ")
|
println("## Random config: ")
|
||||||
println(" - act: ", act)
|
println(" - act: ", act)
|
||||||
println(" - plq: ", plq)
|
println(" - plq: ", plq)
|
||||||
println(" - Epl: ", pl_exact)
|
println(" - Epl: ", pl_exact)
|
||||||
println(" - Ecl: ", cl_exact)
|
println(" - Ecl: ", cl_exact)
|
||||||
|
|
||||||
updt_or_wilson!(U, gp, lp)
|
Ucp = copy(U)
|
||||||
|
|
||||||
|
eo = evenodd(lp)
|
||||||
|
for i in 1:2000
|
||||||
|
updt_or_wilson!(U, gp, eo, lp)
|
||||||
|
end
|
||||||
|
|
||||||
|
plq = plaquette(U, lp, gp, ymws)
|
||||||
|
act = gauge_action(U, lp, gp, ymws)
|
||||||
|
pl_exact = Eoft_plaq(U, gp, lp, ymws)
|
||||||
|
cl_exact = Eoft_clover(U, gp, lp, ymws)
|
||||||
|
println("## After OR [eo]: ")
|
||||||
|
println(" - act: ", act)
|
||||||
|
println(" - plq: ", plq)
|
||||||
|
println(" - Epl: ", pl_exact)
|
||||||
|
println(" - Ecl: ", cl_exact)
|
||||||
|
|
||||||
|
U .= Ucp
|
||||||
|
for i in 1:2000
|
||||||
|
updt_or_wilson!(U, gp, lp)
|
||||||
|
end
|
||||||
|
|
||||||
plq = plaquette(U, lp, gp, ymws)
|
plq = plaquette(U, lp, gp, ymws)
|
||||||
act = gauge_action(U, lp, gp, ymws)
|
act = gauge_action(U, lp, gp, ymws)
|
||||||
|
|
@ -39,3 +62,5 @@ println(" - act: ", act)
|
||||||
println(" - plq: ", plq)
|
println(" - plq: ", plq)
|
||||||
println(" - Epl: ", pl_exact)
|
println(" - Epl: ", pl_exact)
|
||||||
println(" - Ecl: ", cl_exact)
|
println(" - Ecl: ", cl_exact)
|
||||||
|
|
||||||
|
print_timer(linechars = :ascii)
|
||||||
|
|
|
||||||
Loading…
Add table
Add a link
Reference in a new issue