mirror of
https://igit.ific.uv.es/alramos/latticegpu.jl.git
synced 2025-05-14 19:23:42 +02:00
Some cleanup
This commit is contained in:
Alberto Ramos
parent
bee24c52a8
commit
e18183d517
5 changed files with 61 additions and 66 deletions
|
|
@ -31,9 +31,9 @@ Base.one(::Type{SU2{T}}) where T <: AbstractFloat = SU2{T}(one(T),zero(T))
|
|||
Random.rand(rng::AbstractRNG, ::Random.SamplerType{SU2{T}}) where T <: AbstractFloat = exp(SU2alg{T}(randn(rng,T),randn(rng,T),randn(rng,T)))
|
||||
|
||||
"""
|
||||
function normalize(a::SU2)
|
||||
function normalize(a::T) where {T <: Group}
|
||||
|
||||
Return a normalized element of `SU(2)`
|
||||
Return a normalized element of the group.
|
||||
"""
|
||||
function normalize(a::SU2{T}) where T <: AbstractFloat
|
||||
dr = sqrt(abs2(a.t1) + abs2(a.t2))
|
||||
|
|
@ -80,7 +80,7 @@ function isgroup(a::SU2{T}) where T <: AbstractFloat
|
|||
end
|
||||
|
||||
"""
|
||||
function Base.exp(a::SU2alg, t::Number=1)
|
||||
function Base.exp(a::T, t::Number=1) where {T <: Algebra}
|
||||
|
||||
Computes `exp(a)`
|
||||
"""
|
||||
|
|
@ -120,7 +120,7 @@ end
|
|||
|
||||
|
||||
"""
|
||||
function expm(g::SU2, a::SU2alg)
|
||||
function expm(g::G, a::A) where {G <: Algebra, A <: Algebra}
|
||||
|
||||
Computes `exp(a)*g`
|
||||
|
||||
|
|
|
|||
22
src/YM/YM.jl
22
src/YM/YM.jl
|
|
@ -16,12 +16,21 @@ using CUDA, Random, StructArrays
|
|||
using ..Space
|
||||
using ..Groups
|
||||
|
||||
import Base.show
|
||||
|
||||
struct GaugeParm{T}
|
||||
beta::T
|
||||
cG::NTuple{2,T}
|
||||
ng::Int64
|
||||
end
|
||||
export GaugeParm
|
||||
function Base.show(io::IO, gp::GaugeParm)
|
||||
|
||||
println(io, "beta: ", gp.beta)
|
||||
println(io, "Ngauge: ", gp.beta)
|
||||
|
||||
return nothing
|
||||
end
|
||||
|
||||
struct YMworkspace{T}
|
||||
GRP
|
||||
|
|
@ -70,6 +79,19 @@ struct YMworkspace{T}
|
|||
end
|
||||
end
|
||||
export YMworkspace
|
||||
function Base.show(io::IO, ymws::YMworkspace)
|
||||
|
||||
println(io, "Workspace for Group: ", ymws.GRP)
|
||||
println(io, " Algebra: ", ymws.ALG)
|
||||
println(io, "Precision: ", ymws.PRC)
|
||||
if ymws.fpln == nothing
|
||||
println(io, " - Running in memory efficient mode")
|
||||
else
|
||||
println(io, " - Running in computing efficient mode")
|
||||
end
|
||||
return nothing
|
||||
end
|
||||
|
||||
|
||||
include("YMfields.jl")
|
||||
export field, field_pln, randomize!, zero!, norm2
|
||||
|
|
|
|||
|
|
@ -9,29 +9,31 @@
|
|||
### created: Mon Jul 12 18:31:19 2021
|
||||
###
|
||||
|
||||
function krnl_plaq!(plx, U, ipl, lp::SpaceParm)
|
||||
|
||||
id1, id2 = lp.plidx[ipl]
|
||||
b, r = CUDA.threadIdx().x, CUDA.blockIdx().x
|
||||
bu1, ru1 = up((b, r), id1, lp)
|
||||
bu2, ru2 = up((b, r), id2, lp)
|
||||
|
||||
@inbounds plx[b, r] = tr(U[b,id1,r]*U[bu1,id2,ru1] / (U[b,id2,r]*U[bu2,id1,ru2]))
|
||||
|
||||
return nothing
|
||||
end
|
||||
|
||||
function krnl_plaq!(plx, U, lp::SpaceParm)
|
||||
function krnl_plaq!(plx, U, lp::SpaceParm{N,M,D}) where {T,N,M,D}
|
||||
|
||||
b, r = CUDA.threadIdx().x, CUDA.blockIdx().x
|
||||
|
||||
Ush = @cuStaticSharedMem(T, (D,2))
|
||||
|
||||
plx[b,r] = zero(plx[b,r])
|
||||
@inbounds for ipl in 1:lp.npls
|
||||
id1, id2 = lp.plidx[ipl]
|
||||
for id1 in 1:N-1
|
||||
if ru2 == r
|
||||
gt2 = Ush[bu2,1]
|
||||
else
|
||||
gt2 = U[bu2,id1,ru2]
|
||||
end
|
||||
for id2 = id1+1:N
|
||||
if ru1 == r
|
||||
gt1 = Ush[bu1,2]
|
||||
else
|
||||
gt1 = U[bu1,id2,ru1]
|
||||
end
|
||||
sync_threads()
|
||||
|
||||
bu1, ru1 = up((b, r), id1, lp)
|
||||
bu2, ru2 = up((b, r), id2, lp)
|
||||
bu1, ru1 = up((b, r), id1, lp)
|
||||
bu2, ru2 = up((b, r), id2, lp)
|
||||
|
||||
plx[b,r] += tr(U[b,id1,r]*U[bu1,id2,ru1] / (U[b,id2,r]*U[bu2,id1,ru2]))
|
||||
plx[b,r] += tr(Ush[b,1]*gt1 / (Ush[b,2]*gt2))
|
||||
end
|
||||
|
||||
return nothing
|
||||
|
|
@ -139,6 +141,12 @@ function krnl_add_force_plns!(frc::AbstractArray{T}, fpl, lp::SpaceParm{N,M,D})
|
|||
return nothing
|
||||
end
|
||||
|
||||
"""
|
||||
function force_wilson(ymws::YMworkspace, U, lp::SpaceParm)
|
||||
|
||||
Computes the force deriving from the Wilson plaquette action, without
|
||||
the prefactor 1/g0^2, and assign it to the workspace force `ymws.frc1`
|
||||
"""
|
||||
function force_wilson(ymws::YMworkspace, U, lp::SpaceParm)
|
||||
|
||||
if ymws.fpln == nothing
|
||||
|
|
|
|||
|
|
@ -9,8 +9,6 @@
|
|||
### created: Thu Jul 15 15:16:47 2021
|
||||
###
|
||||
|
||||
un(t) = t <: Union{Group, Complex}
|
||||
|
||||
function field(::Type{T}, lp::SpaceParm) where {T}
|
||||
|
||||
sz = lp.bsz, lp.ndim, lp.rsz
|
||||
|
|
|
|||
|
|
@ -9,7 +9,13 @@
|
|||
### created: Thu Jul 15 11:27:28 2021
|
||||
###
|
||||
|
||||
function gauge_action(U, lp::SpaceParm, gp::GaugeParm{T}, ymws::YMworkspace{T}) where T <: AbstractFloat
|
||||
"""
|
||||
|
||||
function gauge_action(U, lp::SpaceParm, gp::GaugeParm, ymws::YMworkspace)
|
||||
|
||||
Returns the value of the Wilson plaquette action for the configuration U.
|
||||
"""
|
||||
function gauge_action(U, lp::SpaceParm, gp::GaugeParm{T}, ymws::YMworkspace{T}) where T <: A
|
||||
|
||||
CUDA.@sync begin
|
||||
CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_plaq!(ymws.cm, U, lp)
|
||||
|
|
@ -64,24 +70,6 @@ function HMC!(U, eps, ns, lp::SpaceParm, gp::GaugeParm, ymws::YMworkspace; noacc
|
|||
return dh, acc
|
||||
end
|
||||
|
||||
function krnl_updt!(mom::AbstractArray{TF}, frc, eps1, U::AbstractArray{TU}, eps2, lp::SpaceParm{N,M,D}) where {TU,TF, N,M,D}
|
||||
|
||||
b, r = CUDA.threadIdx().x, CUDA.blockIdx().x
|
||||
|
||||
Ush = @cuStaticSharedMem(TU, D)
|
||||
Fsh = @cuStaticSharedMem(TF, D)
|
||||
|
||||
@inbounds for id in 1:lp.ndim
|
||||
Ush[b] = U[b,id,r]
|
||||
Fsh[b] = frc[b,id,r]
|
||||
|
||||
mom[b,id,r] = mom[b,id,r] + eps1 * Fsh[b]
|
||||
U[b,id,r] = expm(Ush[b], mom[b,id,r], eps2)
|
||||
end
|
||||
|
||||
return nothing
|
||||
end
|
||||
|
||||
function OMF4!(mom, U, eps, ns, lp::SpaceParm, gp::GaugeParm{T}, ymws::YMworkspace{T}) where T <: AbstractFloat
|
||||
|
||||
r1::T = 0.08398315262876693
|
||||
|
|
@ -94,46 +82,25 @@ function OMF4!(mom, U, eps, ns, lp::SpaceParm, gp::GaugeParm{T}, ymws::YMworkspa
|
|||
ee = eps*gp.beta/gp.ng
|
||||
force_wilson(ymws, U, lp)
|
||||
for i in 1:ns
|
||||
# STEP 1
|
||||
# CUDA.@sync begin
|
||||
# CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_updt!(ymws.mom, ymws.frc1, r1*ee, U, eps*r2, lp)
|
||||
# end
|
||||
mom .= mom .+ (r1*ee) .* ymws.frc1
|
||||
U .= expm.(U, mom, eps*r2)
|
||||
|
||||
# STEP 2
|
||||
force_wilson(ymws, U, lp)
|
||||
# CUDA.@sync begin
|
||||
# CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_updt!(ymws.mom, ymws.frc1, r3*ee, U, eps*r4, lp)
|
||||
# end
|
||||
mom .= mom .+ (r3*ee) .* ymws.frc1
|
||||
U .= expm.(U, mom, eps*r4)
|
||||
|
||||
# STEP 3
|
||||
force_wilson(ymws, U, lp)
|
||||
# CUDA.@sync begin
|
||||
# CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_updt!(ymws.mom, ymws.frc1, r5*ee, U, eps*r6, lp)
|
||||
# end
|
||||
mom .= mom .+ (r5*ee) .* ymws.frc1
|
||||
U .= expm.(U, mom, eps*r6)
|
||||
|
||||
# STEP 4
|
||||
force_wilson(ymws, U, lp)
|
||||
# CUDA.@sync begin
|
||||
# CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_updt!(ymws.mom, ymws.frc1, r5*ee, U, eps*r4, lp)
|
||||
# end
|
||||
mom .= mom .+ (r5*ee) .* ymws.frc1
|
||||
U .= expm.(U, mom, eps*r4)
|
||||
|
||||
# STEP 5
|
||||
force_wilson(ymws, U, lp)
|
||||
# CUDA.@sync begin
|
||||
# CUDA.@cuda threads=lp.bsz blocks=lp.rsz krnl_updt!(ymws.mom, ymws.frc1, r3*ee, U, eps*r2, lp)
|
||||
# end
|
||||
mom .= mom .+ (r3*ee) .* ymws.frc1
|
||||
U .= expm.(U, mom, eps*r2)
|
||||
|
||||
# STEP 6
|
||||
force_wilson(ymws, U, lp)
|
||||
mom .= mom .+ (r1*ee) .* ymws.frc1
|
||||
end
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue