Roadmap

[sostock]

• Simplify tests. Currently, running the tests takes a lot of time. We probably don’t need all of them. (Simplify tests #3)
• Allow knot vectors where the first and last breakpoints are repeated less than k times. ([WIP] Allow more general knot sequences #5)
  • Add a new constructor with mandatory keyword argument:
    • BSplineBasis(k, knots=...) for supplying the knot sequence directly.
    • BSplineBasis(k, breakpoints=...) which creates the knot vector from the breakpoints by repeating the first and last breakpoint k times, equivalent to the old constructor without keyword argument.
  • Indexing a BSplineBasis with an AbstractUnitRange creates a new BSplineBasis with an appropriate knot vector. view(::BSplineBasis, ::AbstractUnitRange) creates an equivalent BSplineBasis whose knot vector is a view. The indices keyword argument can then be removed.
  • Since the BSplineBasis is changed to contain the knot vector instead of the breakpoints, breakpoints(basis) is changed to return unique(knots(basis)). Alternatively, the breakpoints function can be removed completely, since the user can just call unique(knots(basis)) themself.
• Add a StaticBSplineBasis{K, T} type where the order K is a type parameter for potentially better performance.
  • Specialized implementations for small K are possible.
  • Add AbstractBSplineBasis{T} as supertype for BSplineBasis{T} and StaticBSplineBasis{K, T}.

[hyrodium]

Hi! Thank you for this great package!

I'm also developing another package for B-spline: BasicBSpline.jl. With this package, some of the roadmap seem solved.

Allow knot vectors where the first and last breakpoints are repeated less than k times.

As commented here, BasicBSpline.KnotVector can handle any knot vectors.

Add a StaticBSplineBasis{K, T} type where the order K is a type parameter for potentially better performance.

BasicBSpline.jl has a type BasicBSpline.BSplineSpace{p,T} which is for B-spline space with polynomial degree p (== K-1).

With BSplines.jl

julia> using BSplines, BenchmarkTools

julia> basis = BSplineBasis(3, [1,2,3])
4-element BSplineBasis{Vector{Int64}}:
 order: 3
 breakpoints: [1, 2, 3]

julia> length(basis)
4

julia> bsplines(basis,2.1)
3-element OffsetArray(::Vector{Float64}, 2:4) with eltype Float64 with indices 2:4:
 0.4049999999999999
 0.5850000000000001
 0.010000000000000018

julia> @benchmark bsplines(basis,2.1)
BenchmarkTools.Trial: 10000 samples with 967 evaluations.
 Range (min … max):  81.488 ns …   4.399 μs  ┊ GC (min … max): 0.00% … 97.85%
 Time  (median):     85.933 ns               ┊ GC (median):    0.00%
 Time  (mean ± σ):   97.631 ns ± 190.500 ns  ┊ GC (mean ± σ):  9.34% ±  4.68%

   ▁▃▆█▆▇▇▆▅▅▄▄▃▂▂▂▂▁▁▁▁           ▁▁▁▁▁       ▁               ▂
  ▇█████████████████████████▇▇▆▆▇▇███████████▇▇███▇▇▆▆▅▆▇▇█▇▇▄ █
  81.5 ns       Histogram: log(frequency) by time       119 ns <

 Memory estimate: 112 bytes, allocs estimate: 2.

julia> @benchmark bsplines($basis,2.1)
BenchmarkTools.Trial: 10000 samples with 978 evaluations.
 Range (min … max):  66.763 ns …   4.663 μs  ┊ GC (min … max):  0.00% … 98.37%
 Time  (median):     69.476 ns               ┊ GC (median):     0.00%
 Time  (mean ± σ):   82.030 ns ± 209.499 ns  ┊ GC (mean ± σ):  12.50% ±  4.81%

   ▃▅▇█▆▅▄▄▃▃▂▂▂▂▂▁▁                                           ▂
  ▇███████████████████▇█▇▇█▇▆▇▇▇▆▅▅▆▆▇█▇▇▆▇▆▆▆▆▇▇▇▇▇▆▆▅▅▄▆▄▆▇▆ █
  66.8 ns       Histogram: log(frequency) by time       102 ns <

 Memory estimate: 112 bytes, allocs estimate: 2.

With BasicBSplines.jl

julia> using BasicBSpline, BenchmarkTools

julia> k = KnotVector(1,1,1,2,3,3,3)
KnotVector([1.0, 1.0, 1.0, 2.0, 3.0, 3.0, 3.0])

julia> P = BSplineSpace{2}(k)
BSplineSpace{2, Float64}(KnotVector([1.0, 1.0, 1.0, 2.0, 3.0, 3.0, 3.0]))

julia> dim(P)
4

julia> bsplinebasisall(P,2,2.1)
3-element StaticArrays.SVector{3, Float64} with indices SOneTo(3):
 0.4049999999999999
 0.5850000000000001
 0.010000000000000018

julia> @benchmark bsplinebasisall(P,2,2.1)
BenchmarkTools.Trial: 10000 samples with 997 evaluations.
 Range (min … max):  20.290 ns …  1.596 μs  ┊ GC (min … max): 0.00% … 98.44%
 Time  (median):     21.033 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   24.158 ns ± 36.790 ns  ┊ GC (mean ± σ):  3.95% ±  2.59%

  ▅█▇▃             ▁▂▂▃▂         ▁▂                           ▁
  █████▇▆▆▅▅▄▄▅▅▅▇███████▇▆▄▄▄▄▅▄█████▇▆▄▅▄▃▅▅▅▆▄▆▆▇▇▇█▇▇▇▆▄▅ █
  20.3 ns      Histogram: log(frequency) by time      42.9 ns <

 Memory estimate: 32 bytes, allocs estimate: 1.

julia> @benchmark bsplinebasisall($P,2,2.1)
BenchmarkTools.Trial: 10000 samples with 1000 evaluations.
 Range (min … max):  5.710 ns … 16.832 ns  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     5.781 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   5.793 ns ±  0.283 ns  ┊ GC (mean ± σ):  0.00% ± 0.00%

   ▄ ▆▆▂▅█▄ ▂   ▁▃                                           ▂
  ▇█▆██████▁██▄▁██▄▁▄▄▃▃▄▄▃▁▄▅▅▁▄▅▃▄▄▅▄▅▄▅▅▅▄▅▄▅▅▄▅▅▄▅▅▅▅▅▃▃ █
  5.71 ns      Histogram: log(frequency) by time     6.25 ns <

 Memory estimate: 0 bytes, allocs estimate: 0.

If you are interested in BasicBSpline.jl, please check the documentation.

Some of the features in BSplines.jl seem to be missing from BasicBSpline.jl and vice versa. For example, BasicBSpline.jl doesn't have methods for plotting and knotaverages functions.

Merging BSplines.jl and BasicBSpline.jl would be great for other users and developers. Do you have any thoughts on this?

[sostock]

It looks like BasicBSpline.jl does indeed solve the problems I had planned to address at some point. It also does a lot more. I don’t use BSplines.jl anymore and so I haven’t worked on it for a while. (I started to implement arbitrary knot vectors in #5, but that was already a long time ago. I don’t have the time to finish it right now.)

Merging the packages could be nice, but it does not seem straightforward, as the packages do some things quite differently (for example, both packages have a KnotVector type, but they differ strongly in how they behave, even though they have the same purpose). But maybe we can work something out. Do you already have thoughts on how a merged package could look like? Should we start a discussion on how a good API for a merged package could look like? (I’m mostly happy with the API of BSplines.jl, I haven’t looked at BasicBSpline.jl in detail.)

As examples of what is missing in BasicBSpline.jl you mentioned plotting and knotaverages. To me, these seem straightforward to implement. Are there other parts of this package that you would like to see in BasicBSpline.jl? I think your package BasicBSpline.jl can already do most of what I consider essential (and I never intended to implement much more in BSplines.jl). And because it uses static vectors, BasicBSpline.jl is much more performant. To me, it looks like BasicBSpline.jl could mainly profit from a more complete documentation.

[hyrodium]

Thanks for the response!

I started to implement arbitrary knot vectors in #5, but that was already a long time ago. I don’t have the time to finish it right now.

Should we start a discussion on how a good API for a merged package could look like? (I’m mostly happy with the API of BSplines.jl, I haven’t looked at BasicBSpline.jl in detail.)

I’m also happy with the API of BasicBSpline.jl, but I’m not sure how other users think about the API. We both have different API, so I thought it’s a good time to review code with each other. But if you don't have time now, never mind about it.

Are there other parts of this package that you would like to see in BasicBSpline.jl?

Sorry, I haven’t checked the entire feature of BSplines.jl, so I’m not sure. 😅

As examples of what is missing in BasicBSpline.jl you mentioned plotting and knotaverages. To me, these seem straightforward to implement.

To me, it looks like BasicBSpline.jl could mainly profit from a more complete documentation.

Thanks for the suggestion. I’ll deal with them! 😀