Benchmark Julia, Python, R, C

main.jl

# Example adapted from https://github.com/JuliaComputing/JuliaBoxTutorials
# Assumes that Python and R are also installed in the system

# let's compare various implementations of the sum(x) function
x = rand(10^7)
sum(x) # expected result 5×10⁶

#------------------------
# Python implementation
using PyCall

py"""
def py_sum(X):
    s = 0.0
    for x in X:
        s += x
    return s
"""
py_sum = py"py_sum"

py_builtin_sum = pybuiltin("sum")

# test Python implementations
py_sum(x) ≈ sum(x)
py_builtin_sum(x) ≈ sum(x)

#-----------------------
# NumPy implementation

numpy_sum = pyimport("numpy").sum

# test NumPy implementation
numpy_sum(x) ≈ sum(x)

#-------------------
# R implementation
using RCall

R"""
r_sum <- function(X) {
    s = 0.0
    for (i in 1:length(X)) {
        s = s + X[i]
    }
    return(s)
}
"""
r_sum = R"r_sum"

r_builtin_sum = R"sum"

# test R implementations
r_sum(x)[1] ≈ sum(x)
r_builtin_sum(x)[1] ≈ sum(x)

#-----------------------
# Julia implementation

function ju_sum(X)
    s = 0.0
    for x in X
        s += x
    end
    s
end

ju_builtin_sum = sum

function ju_simd_sum(X)
    s = 0.0
    @simd for x in X
        s += x
    end
    s
end

# test Julia implementations
ju_sum(x) ≈ sum(x)
ju_simd_sum(x) ≈ sum(x)

#-------------------
# C implementation
using Libdl

Ccode = """
#include <stddef.h>
double c_sum(size_t n, double *X) {
    double s = 0.0;
    for (size_t i = 0; i < n; ++i) {
        s += X[i];
    }
    return s;
}
"""

# compile to a shared library by piping Ccode to gcc
# (works only if you have gcc installed)
const Clib = tempname() # make a temporary file
open(`gcc -fPIC -O3 -shared -xc -o $(Clib * "." * Libdl.dlext) -`, "w") do f
    print(f, Ccode)
end

# the same as above but with a -ffast-math flag added
const Clib_fastmath = tempname() # make a temporary file
open(`gcc -fPIC -O3 -shared -ffast-math -xc -o $(Clib_fastmath * "." * Libdl.dlext) -`, "w") do f
    print(f, Ccode)
end

# define Julia functions that call the C functions
c_sum(X) = ccall(("c_sum", Clib), Float64, (Csize_t, Ptr{Float64}), length(X), X)

c_fastmath_sum(X) = ccall(("c_sum", Clib_fastmath), Float64, (Csize_t, Ptr{Float64}), length(X), X)

# test C implementations
c_sum(x) ≈ sum(x)
c_fastmath_sum(x) ≈ sum(x)

#---------------
# Benchmarking
using BenchmarkTools

# available implementations
impls = [c_sum, c_fastmath_sum,
         r_sum, r_builtin_sum,
         py_sum, py_builtin_sum, numpy_sum,
         ju_sum, ju_builtin_sum, ju_simd_sum]

# corresponding names
label = ["C", "C -fast-math",
         "R", "R built-in",
         "Python", "Python built-in", "NumPy",
         "Julia", "Julia built-in", "Julia SIMD"]

# loop over implementations
times = map(impls) do impl
    # call function multiple times
    bench = @benchmark ($impl)($x)

    # minimum time in milliseconds
    minimum(bench.times) / 1e6
end

#-----------
# Plotting
using UnicodePlots

# sort times in increasing order
idx   = sortperm(times)
label = label[idx]
times = times[idx]

# times for all implementations
barplot(label, times, title="Time [milliseconds]")

# speedup for all implementations
barplot(label, round.(Int, times[end] ./ times),
        title="Speedup relative to slowest ($(label[end]))")