Skip to content

Conversation

cgarling
Copy link
Contributor

@cgarling cgarling commented Aug 10, 2025

Checklist

  • Appropriate tests were added
  • Any code changes were done in a way that does not break public API
  • All documentation related to code changes were updated
  • The new code follows the
    contributor guidelines, in particular the SciML Style Guide and
    COLPRAC.
  • Any new documentation only uses public API

Additional context

This hopes to address the performance issues of multi-valued interpolation discussed in #431 by using dependent variable arrays with elements of StaticArrays.SArray to remove allocations and improve performance. On master this still allocates when evaluating the interpolants, hurting performance.

For the example of CubicSpline, I isolated the allocations to come from the z array. By converting the elements of z to match eltype(u), we are able to remove allocations when evaluating the interpolant and greatly improve the performance. Example usage and performance comparison is given below.

All tests pass for me locally on Julia 1.11.6. If this seems like a reasonable solution, I can look into modifying the other relevant methods as well.

using StaticArrays: SVector, SArray
using DataInterpolations: CubicSpline
using BenchmarkTools: @benchmark
using  Test: @test

x = 1:1_000
a,b,c = [rand(length(x)) for i in 1:3]
y = [SVector{3, Float64}(a[i], b[i], c[i]) for i in eachindex(a,b,c)]

cubic = CubicSpline(y, x)
@benchmark $cubic($2.5)
# On master, ↑ is 133.750 ns with Memory estimate: 640 bytes, allocs estimate: 16.
# On PR: 21.717 ns, 0 allocations

y2 = [SArray{Tuple{1,2,3}, Float64, 3, 6}(rand(2,3)) for i in eachindex(a,b,c)]
cubic2 = CubicSpline(y2, x)
@benchmark $cubic2($2.5)
# On master, ↑ is 248.704 ns with Memory estimate: 1.00 KiB, allocs estimate: 16.
# On PR: 25.965 ns ns, 0 allocations

# Vector of vectors
y3 = [[a[i], b[i], c[i]] for i in eachindex(a,b,c)]
cubic3 = CubicSpline(y3, x)
@benchmark $cubic3($2.5)
# On master, ↑ is 243.811 ns with Memory estimate: 1.17 KiB, allocs estimate: 30.
# On PR: 229.563 ns with Memory estimate: 1.17 KiB, allocs estimate: 30.

# For Matrix y
y4 = [a b c]'
x4 = 1:size(y4, 2)
cubic4 = CubicSpline(y4, x4)
@benchmark $cubic4($2.5)
# On master, ↑ is 307.918 ns with Memory estimate: 1.64 KiB, allocs estimate: 42.
@test cubic(2.5) == cubic4(2.5) # Gives same result

@ChrisRackauckas
Copy link
Member

It looks reasonable enough to me. Though what does it do to the arrays of arrays case? allocate the views?

@cgarling
Copy link
Contributor Author

cgarling commented Aug 10, 2025

I think so, results for y::Vector{Vector} (named y3) added above.

If you meant the CubicSpline(u::AbstractArray{T, N}, t) method that can be used for multi-valued interpolation, then also yes -- y4 above shows a matrix input y4 and vector x4 for multi-valued interpolation, but it is much slower than the implementation with y::Vector{SVector}.

Covers vector of vector input and vector of matrices input
Can support `Vector{Vector}` and `Vector{Matix}` input as well; tests added to confirm behavior.
Need to remove value `push!`ed onto `u` prior to new tests.
@cgarling
Copy link
Contributor Author

With this PR, the following methods work with u::AbstractVector{<:AbstractVector}, and seem to have good performance and no allocations when u::AbstractVector{<:StaticArrays.SVector}.

  • CubicSpline
  • QuadraticInterpolation
  • SmoothedConstantInterpolation
  • ConstantInterpolation
  • QuadraticSpline
  • LinearInterpolation
  • CubicHermiteSpline
  • QuinticHermiteSpline

Not all of these methods have tests confirming this behavior. I had to relax the type signatures of the _interpolate methods for CubicHermiteSpline and QuinticHermiteSpline to enable this behavior, so I added tests confirming that the alterations result in successful operation. I could look at adding similar extra tests to the other interpolation types above, if you'd like to have this behavior covered by the tests, but that should probably be under a different PR.

The remaining methods that don't seem to work with u::AbstractVector{<:AbstractVector} are

  • PCHIPInterpolation
  • AkimaInterpolation
  • LagrangeInterpolation
  • BSplineInterpolation
  • BSplineApprox

I am interested in trying to get PCHIPInterpolation and AkimaInterpolation working, so I may look into that in another PR. For now I think this PR is ready for review.

@cgarling cgarling marked this pull request as ready for review August 11, 2025 12:30
@ChrisRackauckas
Copy link
Member

and seem to have good performance and no allocations when u::AbstractVector{<:StaticArrays.SVector}.

We should test this via AllocCheck.jl so that it doesn't regress.

@cgarling
Copy link
Contributor Author

I have pushed an example for what tests checking StaticArray behavior may look like. If this looks good I can make similar tests for the other interpolations

@ChrisRackauckas ChrisRackauckas merged commit befa6ac into SciML:master Sep 3, 2025
17 of 19 checks passed
@cgarling cgarling deleted the staticarrays2 branch September 3, 2025 13:18
@cgarling cgarling mentioned this pull request Sep 3, 2025
5 tasks
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants