Performance Annotations

Julia provides several macros and annotations that can help improve performance by giving the compiler additional information or relaxing certain constraints.

Fast Math

We can allow the compiler to include floating point optimisations that are correct for real numbers, but may lead to differences for IEEE encoded floats. These may change numerical results and accuracy, but may improve the performance of your code. This is enabled via the @fastmath macro in Julia:

function nofastmath_example!(y, x)
    @inbounds for i in eachindex(y, x)
        y[i] = sin(sqrt(x[i] * x[i] + 1.0) / cos(x[i] + 0.1))
    end
    nothing
end

function fastmath_example!(y, x)
    @inbounds @fastmath for i in eachindex(y, x)
        y[i] = sin(sqrt(x[i] * x[i] + 1.0) / cos(x[i] + 0.1))
    end
    nothing
end

We can benchmark these two algorithms:

import BenchmarkTools: @benchmark, @btime, @belapsed

x = rand(1024) .+ 1.0; y = similar(x);
display(@benchmark nofastmath_example!($y, $x))

BenchmarkTools.Trial: 10000 samples with 1 evaluation per sample.
 Range (min … max):  32.150 μs … 251.150 μs  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     34.330 μs               ┊ GC (median):    0.00%
 Time  (mean ± σ):   35.192 μs ±   5.010 μs  ┊ GC (mean ± σ):  0.00% ± 0.00%

  ▁▂▁▄███▅▆▇▃▄        ▁▁▁▁                                     ▂
  ███████████████▇▇▇▇███████▇▇▇▇▆▆▆▆▅▆▅▄▄▄▄▄▁▁▄▄█▄▁▁▁▁▆▅▄▆▅▃▃▇ █
  32.2 μs       Histogram: log(frequency) by time        53 μs <

 Memory estimate: 0 bytes, allocs estimate: 0.

display(@benchmark fastmath_example!($y, $x))

BenchmarkTools.Trial: 10000 samples with 1 evaluation per sample.
 Range (min … max):  32.190 μs … 64.880 μs  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     33.800 μs              ┊ GC (median):    0.00%
 Time  (mean ± σ):   34.184 μs ±  1.371 μs  ┊ GC (mean ± σ):  0.00% ± 0.00%

         ▂▆██▇▆▄▄▆▇▃   ▂                                      ▂
  ██▆▇▅▁▇████████████▆███▆▇▇▆▇▅▆▇▅▆▆▄▃▅▃▅▅▆▅▇▅▅▇▇▇▇█▇███▇▇▇▇▇ █
  32.2 μs      Histogram: log(frequency) by time      40.8 μs <

 Memory estimate: 0 bytes, allocs estimate: 0.

On this machine the two versions come out essentially identical, with @fastmath giving no measurable improvement at all. That is not unusual: the dominant cost in this loop is the sin and cos calls, and relaxing IEEE semantics has very little to offer on those. It is still useful to keep in mind, as on some systems and for some expressions it lets the compiler use faster floating-point instructions and reorder computations.

Using the “fast” math operations violates strict IEEE semantics. The compiler is allowed to assume that no argument is a NaN or an Inf, to reassociate sums and products, and to replace a division by a reciprocal multiply. This means @fastmath can change the answer your program produces, and can silently break code that relies on NaN propagation. For this reason it is often avoided in scientific applications, and is an opt-in performance enhancement that you should benchmark and validate rather than sprinkle over a codebase.

Bounds Checking

We have already seen the use of the @inbounds macro throughout this course. This is one of the easiest optimisations to make, as long as you are confident that you are accessing memory in a correct way. Turning off bounds checking and accessing incorrect areas of memory may lead to undefined behaviour, memory corruption and crashes. This can be mitigated by proper use of methods like eachindex or axes.

function with_bounds_check(x)
    s = zero(eltype(x))
    for i in 1:length(x)
        s += x[i]
    end
    s
end

function without_bounds_check(x)
    s = zero(eltype(x))
    @inbounds for i in 1:length(x)
        s += x[i]
    end
    s
end

x = rand(1000);
display(@benchmark with_bounds_check($x))

BenchmarkTools.Trial: 10000 samples with 147 evaluations per sample.
 Range (min … max):  694.966 ns …  3.874 μs  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     707.551 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   720.540 ns ± 73.636 ns  ┊ GC (mean ± σ):  0.00% ± 0.00%

  ▄ █                                                           
  █▇█▆█▃▃▃▃▃▃▂▂▂▂▂▂▂▂▂▂▁▂▂▁▁▂▁▁▁▁▁▂▁▁▁▁▂▁▁▂▁▁▁▁▂▂▁▁▁▂▁▁▁▁▁▁▁▂▂ ▂
  695 ns          Histogram: frequency by time         1.06 μs <

 Memory estimate: 0 bytes, allocs estimate: 0.

display(@benchmark without_bounds_check($x))

BenchmarkTools.Trial: 10000 samples with 151 evaluations per sample.
 Range (min … max):  683.179 ns …  1.894 μs  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     708.874 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   721.676 ns ± 66.763 ns  ┊ GC (mean ± σ):  0.00% ± 0.00%

    ▇█▄▃                                                        
  █▅█████▄▄▃▃▃▃▃▂▂▂▃▂▂▂▂▂▂▂▂▂▂▂▂▂▂▂▂▁▁▂▂▁▂▁▁▂▂▂▂▂▂▂▂▁▁▁▁▁▁▂▂▁▂ ▃
  683 ns          Histogram: frequency by time         1.06 μs <

 Memory estimate: 0 bytes, allocs estimate: 0.

In this particular case the two are indistinguishable: the loop is simple enough that the compiler can hoist the single bounds check out of the loop, so there is nothing left to remove. The improvement from removing bounds checking can still be significant in tight loops where the compiler cannot prove the accesses are safe. This optimisation was more necessary in older versions of Julia, but is not always required, especially when combined with semantics such as eachindex or axes. Note also that @inbounds applied to an index that is actually out of range is undefined behaviour — Julia will happily read or write memory that does not belong to the array, and the result may be a wrong answer rather than a crash.

Inlining

If a function is short enough, one can encourage the compiler to always inline the function so that at runtime, one does not have to pay the cost of calling another function. This can be done in Julia using the @inline macro:

@inline function fast_add(a, b)
    a + b
end

function use_fast_add(x)
    s = zero(eltype(x))
    for val in x
        s = fast_add(s, val)
    end
    s
end

This makes a strong suggestion to the compiler that wherever fast_add is used, it should be directly inlined, instead of inserting a function call at the call site. We should note that the compiler contains a heuristic for whether or not to inline your code for you, and using the @inline macro is only a suggestion to the compiler to try and encourage inlining. Most of the time this is not necessary, but might be useful in performance critical code.

Constant Propagation

A compiler can remove code entirely if all the constants involved are known at compile time, evaluating the result once during compilation instead of on every call. Take the following example, using the @code_typed macro from the InteractiveUtils standard library (loaded automatically in the REPL, but needing an explicit import inside a script):

using InteractiveUtils
compiled_fn() = sum(1:1000);
display(@code_typed compiled_fn())

CodeInfo(
1 ─     return 500500
) => Int64

One can see that the actual code simply returns the constant value which was calculated. The code never actually performs the sum at runtime, since it can be precomputed at compile time.

If you can give the compiler information about constants during compile time, it can propagate that information forwards to avoid costly computations down the line. There is a special type in Julia called Val which allows us to insert data into the type information:

function dynamics_rule_with_val(u, ::Val{N}) where {N}
    unit = one(typeof(u))
    mask = ~(~zero(typeof(u)) << N)
    u_left = (u << 1) | ((u & (unit << (N-1))) >> (N-1))
    u_right = (u >> 1) | ((u & unit) << (N-1))
    return (xor(xor(u_left, u), u_right) & mask)
end

function dynamics_rule_no_val(u, N)
    unit = one(typeof(u))
    mask = ~(~zero(typeof(u)) << N)
    u_left = (u << 1) | ((u & (unit << (N-1))) >> (N-1))
    u_right = (u >> 1) | ((u & unit) << (N-1))
    return (xor(xor(u_left, u), u_right) & mask)
end

u = rand(Int)
N = 32
display(@benchmark dynamics_rule_no_val($u, $N))

BenchmarkTools.Trial: 10000 samples with 1000 evaluations per sample.
 Range (min … max):  3.810 ns … 14.660 ns  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     3.860 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   3.925 ns ±  0.494 ns  ┊ GC (mean ± σ):  0.00% ± 0.00%

   █                                                          
  ▃█▆▅▅▂▂▂▂▂▂▁▂▂▂▁▂▂▂▂▂▂▂▂▂▂▂▂▂▁▂▂▂▂▂▂▂▁▂▂▂▂▂▂▂▂▁▂▂▂▂▂▂▂▁▂▂▂ ▂
  3.81 ns        Histogram: frequency by time        5.29 ns <

 Memory estimate: 0 bytes, allocs estimate: 0.

display(@benchmark dynamics_rule_with_val($u, $(Val(N))))

BenchmarkTools.Trial: 10000 samples with 1000 evaluations per sample.
 Range (min … max):  2.350 ns … 48.270 ns  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     2.420 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   2.456 ns ±  0.716 ns  ┊ GC (mean ± σ):  0.00% ± 0.00%

  ▃   ▅█▃  ▅▁                                                ▁
  █▃▁▁███▄▆██▁▄▁▁▁▄▄▃▁▄▅▄▃▃▄▄▄▄▃▄▄▃▅▁▄▄▄▃▄▃▅▅▃▅▅▄▅▄▄▅▅▅▅▅▆▆▆ █
  2.35 ns      Histogram: log(frequency) by time     3.04 ns <

 Memory estimate: 0 bytes, allocs estimate: 0.

Here the Val version is noticeably faster — it removes roughly a third of the runtime — because N is part of the type, so mask and the shift amounts are folded into constants at compile time rather than being recomputed on every call. Val can also be used to try and make your code as generic as possible, while still retaining as much performance as possible. This technique should only be used if the data contained within the Val type only has a few possible values, and likely won’t change during the execution of your code.

Note: When introducing Val into your code, you likely have to propagate this changed call signature up the call stack, which might involve significant changes to your codebase. This is because constructing a Val type is usually not type-safe (not predictable by the compiler), and so you want to do this as far up the call stack as possible to avoid the poor performance of type unstable code. In these cases, it is sometimes better to just pass in the datatype you want to use.