Comparison with LinearInterpolators

[JeffFessler]

For the comparison with LinearInterpolators.jl,
I haven't got cubic B-spline to work yet due to a Compat issue.

But here is the linear B-spline comparison:

using LinearInterpolators: SparseInterpolator, LinearSpline
# using InterpolationKernels: BSpline # todo later
using LinearAlgebra: mul!
using Plots
using Random: seed!
default(markerstrokecolor = :auto, label="")
seed!(0)

T = Float64
#T = Float32
N = 100
xs = range(T(0), T(10), N) # source grid
ys = sin.(xs) # initial y data

M = 1401
xf = range(T(-2), T(12), M) # fine grid, just for plotting
yf = sin.(xf)

Q = 1000
xq = sort(14 * rand(T, Q) .- 2) # query points for interpolation
yq = sin.(xq)

p1 = plot()
plot!(xf, yf, label="fine")
scatter!(xs, ys, label="data", marker=:x)
scatter!(xq, yq, label="query", marker=:o)

ker = LinearSpline(T)
# ker = BSpline{4,T} # cubic !
interp1 = SparseInterpolator(ker, xq, xs)

out1 = similar(yq)
mul!(out1, interp1, ys)

plot!(p1, xq, out1, label="interp1")


using FastInterpolations: linear_interp!
using FastInterpolations: ClampExtrap
using FastInterpolations: ExtendExtrap

out2 = similar(yq)

#extrap = ExtendExtrap()
extrap = ClampExtrap()
linear_interp!(out2, xs, ys, xq; extrap)

plot!(p1, xq, out2, label="interp2")
gui()

@assert out1 ≈ out2


using BenchmarkTools
@btime mul!($out1, $interp1, $ys); # 2.2 us - 0 alloc
@btime linear_interp!($out2, $xs, $ys, $xq; extrap); # 2.8 us, 1 alloc 64 byte

#=
with T = Float32, I get 
1.7 μs (0 allocations: 0 bytes) - mul!
2.9 μs (1 allocation: 48 bytes) - linear_interp!
=#
nothing

For my work Float32 suffices and curiously, FastInterpolations.jl was slighly slower for 32bit than 64bit, whereas the other package is faster. I had expected both to be faster - less memory to move around.

versioninfo()

Julia Version 1.12.5
Commit 5fe89b8ddc1 (2026-02-09 16:05 UTC)
Build Info:
  Official https://julialang.org release
Platform Info:
  OS: macOS (x86_64-apple-darwin24.0.0)
  CPU: 8 × Intel(R) Core(TM) i7-7700K CPU @ 4.20GHz
  WORD_SIZE: 64
  LLVM: libLLVM-18.1.7 (ORCJIT, skylake)
  GC: Built with stock GC
Threads: 8 default, 1 interactive, 8 GC (on 8 virtual cores)
Environment:
  JULIA_NUM_THREADS = 8

[mgyoo86]

Hi, thanks for the detailed comparison.

I think what you observed comes from benchmarking artifacts caused by how the extrap keyword argument was passed, which introduced extra allocation (1 alloc) and so the noise. It would be better to interpolate the keyword argument with $ as well, just like the other benchmark inputs:

# Before
# This allocates because `extrap` is captured
@btime linear_interp!($out2, $xs, $ys, $xq; extrap);

# Recommended
# This benchmarks the call more accurately with zero allocation
@btime linear_interp!($out2, $xs, $ys, $xq; extrap=$extrap);

On my local machine, I found that FastInterpolations.jl is slightly faster than LinearInterpolators.jl for both Float32 and Float64.

That said, linear interpolation itself is such a simple operation that there may not be much room for further optimization.
In practice, the bottleneck is often not the floating-point arithmetic itself, but other overheads such as memory access, interval search, and function-call overhead. So while Float32 can in principle be a bit faster than Float64, in this kind of interpolation benchmark the difference may be very small or even within the noise.

It also seems that LinearInterpolators.jl may not support vector grids properly. When I tested that case, I ran into a stack overflow.

Please try rerunning the benchmark with the corrected setup, and let me know if you still see anything unusual.

[JeffFessler]

Oops, sorry I didn't realize that even a scalar keyword argument should be interpolated. Your modification eliminated the alloc, and then I got machine-dependent timing results.
On a 2020 Intel iMac, I got 1.6us for LinearInterpolators and 2.5us for FastInterpolations.
On a 2023 Apple M3Max I got 1.0us for both. In both cases the 32 and 64-bit results took the same amount of time. (Maybe it would matter more for 2D interpolation of bigger arrays.). Clearly the newer machine's results are more relevant, so this alone is enough to convince me to switch to FastInterpolations for future work. Still hoping to eventually compare the cubic B-spline case for completeness.

(My uses cases all use uniform grids for input data, so I have never looked at any vector grids.)

[mgyoo86]

It is interesting that the timings are strongly machine-dependent than I expected.
Anyway, I'm glad to hear that the results on the newer machine look good.

High performance is the top priority for this package,
so if you notice any similar regressions in the future, please let me know :)

It would also be nice to see a cubic-spline comparison at some point, if feasible.
That could also make a useful addition to the benchmark graph in README as you suggested before.

[mgyoo86]

@JeffFessler
I've revisited this issue after the v0.4.5 release.

The machine-dependent performance issue you observed (on Intel iMac)
turned out to be the cost of TwicePrecision arithmetic, which is much more expensive on Intel architecture than Apple M-series architecture (see Issue #85 and PR #87 for details). That has now been addressed, and I also reduced a few other small overheads, so I ran the benchmark again.

Here are the current results on my side:

CPU	LinearInterpolators	FastInterpolations.jl	Result
Apple M1	1.5 μs	0.971 μs	~50% faster
Intel Xeon	2.0 μs	1.78 μs	~12% faster

From this quick benchmark, FastInterpolations.jl seems to be consistently faster than LinearInterpolators.jl for the linear interpolation.

[JeffFessler]

Oh wow, thanks so much for telling me about that. I had no idea that the TwicePrecision would have that effect. It's surprising that you had to introduce a new type. Seems like Julia would have already had some way to control the precision.

The LinRange type doesn't use TwicePrecision, which is ironic because I switched from using LinRange to range in all my code some time ago, just for less typing :)

I was able to essentially reproduce your results with my 2020 intel iMac:

Float32: 1.4 us LinearInterpolators | 1.6 us FastInterpolations
Float64: 1.7 us LinearInterpolators | 1.5 us FastInterpolations