pi_avx_transducers.jl

using InteractiveUtils, Pkg

for pk in ["Transducers", "LoopVectorization"]
    if pk ∉ keys(Pkg.project().dependencies)
        Pkg.add(pk)
    end
end

using Transducers, LoopVectorization
function picalc_kernel(r, a, b)
    su = 0.0
    @avx for j in a:b
        su += 4.0/(1.0 + r[j]*r[j])
    end
    su
end

function picalc(numsteps)
    println("Calculating PI using:")
    println("  ", numsteps, " slices")
    w = Threads.nthreads()
    println("  ", w, " threads(s)")
    start = time()
    slice = 1.0/numsteps
    r = ((1:numsteps) .- 0.5)*slice;
    step = numsteps÷w
    @assert numsteps%w == 0

    sum = reduce(+, Map(k -> picalc_kernel(r, k, k+step-1)), 1:step:numsteps, basesize=1)
    mypi = sum * slice
    
    stop = time()
    elapsed = (stop - start)
    println("Obtained value of PI: ", mypi)
    println("Time taken: ", elapsed, " seconds")
    mypi
end

numsteps=1000000000

@info "Warmup"
picalc(numsteps ÷ 100000)

println()

@info "Real shit"
picalc(numsteps)

println()

versioninfo()

println()