Discussions for Gaugefields.jl in English

[cometscome]

Questions, comments, and bug reports related to LatticeQCD.jl in English. If you are unsure about the category or how to use it, please post here.

[phdtutorhub]

Is MPI broken in Gaugefields.jl after it is updated to incorporate CUDA and JACC? Is JACC in the 't Hooft twisted b.c./flux version of Gaugefields.jl still work in progress? If so, any idea when MPI and JACC in the 't Hooft version of Gaugefields.jl be fixed?

[akio-tomiya]

Hi. Thanks for using Guaugefields.jl. I will ask a contributer who is working it.
Can you tell us what kind of error you find?
Sorry for inconvenience.

[phdtutorhub]

Hi, thanks for the quick reply.

I was running the sample HMC with MPI code on https://github.com/akio-tomiya/Gaugefields.jl. MPI was working before until version v0.6.2 and v0.6.3. It is possible that the inclusion of CUDA and JACC might have broken MPI. It could be a minor bug.

When I run the code with JACC and a CUDA backend on my laptop (mpiexecjl -np 2 -ppn 2 julia --project /home/atang/slhmc/test/gaugefields.jl 1 1 1 2 true), I got this error:
[ Info: CUDA backend loaded
[ Info: CUDA backend loaded
Vector{Gaugefields.AbstractGaugefields_module.Gaugefields_4D_nowing_mpi{3}}
ERROR: ERROR: LoadError: LoadError: type Gaugefields_4D_nowing_mpi has no field win_i
Stacktrace:
type Gaugefields_4D_nowing_mpi has no field win_i

When I run the same code on CPU without JACC, I got the same error:
Vector{Gaugefields.AbstractGaugefields_module.Gaugefields_4D_nowing_mpi{3}}
ERROR: LoadError: ERROR: LoadError: type Gaugefields_4D_nowinheg_mpi has no field win_i

But if I run these codes without MPI:
mpiexecjl -np 2 -ppn 2 julia --project /home/atang/slhmc/test/gaugefields.jl 1 1 1 2 false
mpiexecjl -np 2 -ppn 2 julia --project /home/atang/slhmc/test/gaugefields-cpu.jl 1 1 1 2 false
they work. It is obvious that the problem is MPI. More specifically, Gaugefields_4D_nowing_mpi needs to have a field for win_i.

I also tried running MoriKawa san's version of Gaugefields.jl with 't Hooft twisted bounday conditions using his HMC with MPI code snippet on https://github.com/o-morikawa/Gaugefields.jl. I have modified his code snippet in order for it to run. My modification is appended at the end of this message. When I run his code with JACC and MPI (mpiexecjl -np 2 -ppn 2 julia --project /home/atang/slhmc/test/tflux.jl 1 1 1 2 true), I got this error:
[ Info: CUDA backend loaded
[ Info: CUDA backend loaded
ERROR: ERROR: LoadError: LoadError: UndefVarError: thooftFlux_4D_B_at_bndry_nowing_mpi not defined in Gaugefields.Bfield_module

I got the same error when I ran the code on CPU without JACC (mpiexecjl -np 2 -ppn 2 julia --project /home/atang/slhmc/test/tflux-cpu.jl 1 1 1 2 true):
ERROR: LoadError: ERROR: LoadError: UndefVarError: thooftFlux_4D_B_at_bndry_nowing_mpi not defined in Gaugefields.Bfield_module

But when I ran the codes without MPI:
mpiexecjl -np 2 -ppn 2 julia --project /home/atang/slhmc/test/tflux.jl 1 1 1 2 false
mpiexecjl -np 2 -ppn 2 julia --project /home/atang/slhmc/test/tflux-cpu.jl 1 1 1 2 false
they worked. Again the problem is with MPI. It is probably a minor bug.

At last, Initialize_Bfields in Morikawa san's implementation does not accept the accelerator argument so that it cannot be run with JACC. I assume that accelerator is an important feature to include.

Below is my version of HMC with MPI with 't Hooft twisted boundary conditions:

using Random
using Gaugefields
using Wilsonloop
using LinearAlgebra
using MPI
using JACC

JACC.@init_backend

if length(ARGS) < 5
error("USAGE: ","""
mpirun -np 2 exe.jl 1 1 1 2 true
""")
end
const pes = Tuple(parse.(Int64,ARGS[1:4]))
const mpi = parse(Bool,ARGS[5])

function calc_action(gauge_action,U,p)
NC = U[1].NC
Sg = -evaluate_GaugeAction(gauge_action,U)/NC #evaluate_Gauge_action(gauge_action,U) = tr(evaluate_Gaugeaction_untraced(gauge_action,U))
Sp = p*p/2
S = Sp + Sg
return real(S)
end

function MDstep!(gauge_action,U,p,MDsteps,Dim,Uold,temp1,temp2)
Δτ = 1.0/MDsteps
gauss_distribution!(p)
Sold = calc_action(gauge_action,U,p)
substitute_U!(Uold,U)

for itrj=1:MDsteps
    U_update!(U,p,0.5,Δτ,Dim,gauge_action)

    P_update!(U,p,1.0,Δτ,Dim,gauge_action,temp1,temp2)

    U_update!(U,p,0.5,Δτ,Dim,gauge_action)
end
Snew = calc_action(gauge_action,U,p)
if get_myrank(U) == 0
    println("Sold = $Sold, Snew = $Snew")
    println("Snew - Sold = $(Snew-Sold)")
end
ratio = min(1,exp(-Snew+Sold))
r = rand()
if mpi
    r = MPI.bcast(r, 0, MPI.COMM_WORLD)
end
#ratio = min(1,exp(Snew-Sold))
if r > ratio
    substitute_U!(U,Uold)
    return false
else
    return true
end

end

function U_update!(U,p,ϵ,Δτ,Dim,gauge_action)
temps = get_temporary_gaugefields(gauge_action)
temp1 = temps[1]
temp2 = temps[2]
expU = temps[3]
W = temps[4]

for μ=1:Dim
    exptU!(expU,ϵ*Δτ,p[μ],[temp1,temp2])
    mul!(W,expU,U[μ])
    substitute_U!(U[μ],W)

end

end

function P_update!(U,p,ϵ,Δτ,Dim,gauge_action,temp1,temp2) # p -> p +factorUdSdUμ
NC = U[1].NC
temp = temp1
dSdUμ = temp2
factor = -ϵ*Δτ/(NC)

for μ=1:Dim
    calc_dSdUμ!(dSdUμ,gauge_action,μ,U)
    mul!(temp,U[μ],dSdUμ) # U*dSdUμ
    Traceless_antihermitian_add!(p[μ],factor,temp)
end

end

function HMC_test_4D_tHooft(NX,NY,NZ,NT,NC,Flux,β)
Dim = 4
Nwing = 0

flux = Flux

Random.seed!(123)

if mpi
    PEs = pes#(1,1,1,2)
    U = Initialize_Gaugefields(NC,Nwing,NX,NY,NZ,NT,condition = "hot",mpi=true,PEs = PEs,mpiinit = false,accelerator="JACC") 
    B = Initialize_Bfields(NC,flux,Nwing,NX,NY,NZ,NT,condition = "tflux",mpi=true,PEs = PEs,mpiinit = false)
else
    U = Initialize_Gaugefields(NC,Nwing,NX,NY,NZ,NT,condition = "hot",accelerator="JACC")
    B = Initialize_Bfields(NC,flux,Nwing,NX,NY,NZ,NT,condition = "tflux")
end

if get_myrank(U) == 0
    println("Flux : ", flux)
end

if get_myrank(U) == 0
    println(typeof(U))
end


#temps = Temporalfields(U[1], num=10)
#comb, factor = set_comb(U, Dim)

if Dim == 4
    comb = 6 #4*3/2
elseif Dim == 3
    comb = 3
elseif Dim == 2
    comb = 1
else
    error("dimension $Dim is not supported")
end

factor = 1/(comb*U[1].NV*U[1].NC)

temp1 = similar(U[1])
temp2 = similar(U[1])

#@time plaq_t = calculate_Plaquette(U,B,temp1,temp2)*factor
@time plaq_t = calculate_Plaquette(U,temp1,temp2)*factor
if get_myrank(U) == 0
    println("0 plaq_t = $plaq_t")
end
poly = calculate_Polyakov_loop(U,temp1,temp2) 
if get_myrank(U) == 0
    println("0 polyakov loop = $(real(poly)) $(imag(poly))")
end

#gauge_action = GaugeAction(U,B)
gauge_action = GaugeAction(U)
plaqloop = make_loops_fromname("plaquette")
append!(plaqloop,plaqloop')
β = β/2
push!(gauge_action,β,plaqloop)

#show(gauge_action)

p = initialize_TA_Gaugefields(U) #This is a traceless-antihermitian gauge fields. This has NC^2-1 real coefficients. 
Uold = similar(U)
MDsteps = 50
numaccepted = 0

numtrj = 100
for itrj = 1:numtrj
    t = @timed begin
        accepted = MDstep!(gauge_action,U,p,MDsteps,Dim,Uold,temp1,temp2)
    end
    if get_myrank(U) == 0
        println("elapsed time for MDsteps: $(t.time) [s]")
    end
    numaccepted += ifelse(accepted,1,0)

    #plaq_t = calculate_Plaquette(U,B,temps)*factor
    #println("$itrj plaq_t = $plaq_t")
    
    if itrj % 10 == 0
        #plaq_t = calculate_Plaquette(U,B,temp1,temp2)*factor
    plaq_t = calculate_Plaquette(U,temp1,temp2)*factor
        if get_myrank(U) == 0
            println("$itrj plaq_t = $plaq_t")
        end
        poly = calculate_Polyakov_loop(U,temp1,temp2) 
        if get_myrank(U) == 0
            println("$itrj polyakov loop = $(real(poly)) $(imag(poly))")
            println("acceptance ratio ",numaccepted/itrj)
        end
    end
end
return plaq_t,numaccepted/numtrj

end

function main()
β = 5.7
NX = 4
NY = 4
NZ = 4
NT = 4
NC = 3
Flux = [0,0,1,1,0,0]
#HMC_test_4D(NX,NY,NZ,NT,NC,β)
HMC_test_4D_tHooft(NX,NY,NZ,NT,NC,Flux,β)
end
main()

[cometscome]

Hi. Thank you for sharing your code.
I am Yuki Nagai, one of the developers of JuliaQCD. In your code, you try to use MPI AND CUDA.
We did not implement multi-GPU code so Gaugefield with both mpi and CUDA (or CUDA) is not supported yet.
Please use mpi or JACC or CUDA separately.
In the case of the MPI only version, I found a bug so please use v6.8 or above.
Now I am implementing a hybrid parallel code in v0.6. If you use MPILattice flag, you can use both MPI and GPU. But I would like to note that multi-GPU code is slow. So single GPU or MPI with threads are better to use.
A sample code is shown in
https://github.com/akio-tomiya/Gaugefields.jl/blob/master/test/MPIJACCtest/mpitest.jl