3
1
B-splines are defined using a set of "knots", "control points", and degree. For the purposes of this challenge, the control points are not needed.
The knots are simply a list of strictly increasing real values (note that real B-splines may have repeated knots, for the purposes of this challenge assume knots are all unique). The degree is a non-negative integer.
The set of polynomial basis functions are recursively defined as (image credit):
where i
is the i-th
basis function, and p
is the degree of the basis function. u
is the curve parameterization, and u_j
is the j'th
knot.
The goal of this challenge is to evaluate all the basis functions for a curve of degree m
at location u
.
Input
Your program/function should take as input 3 parameters:
- A real value for the location
u
- A non-negative integer for the curve degree
m
- A list of strictly increasing real values for the knots vector
u_j
You may take these parameters in any order desired and may mutate the inputs as long as no additional information is added (for example, you may take the knots vector to be a combination of a real pointer and integer length).
The input may come from any source desired (stdio, function parameter, etc.).
You may assume the input describes a valid B-spline and that u
is in the range of the knots vector.
Output
Your program/function should output a list of real numbers describing what the basis functions evaluate to. That is, you should output for all
0 <= i <= n-1-m
, where n
is the number of knots (For 1-based indexing, then 1 <= i <= n - m
). The output must be accurate to within 1%
or 1000 ULP of the true value, whichever is greater. The output may be to any sink desired (stdout, return value, etc.)
Examples
All examples are formatted as such:
u
m
u_0 u_1 u_2 ...
N_{0,m}(u) N_{1,m}(u) N_{2,m}(u) ...
case 1:
.5
0
0 0.5 1
0 1
case 2:
.5
1
-1 0 1 2
0.5 0.5
case 3:
0.9
1
-1 0 1 2 3
0.1 0.9 0
case 4:
1.1
1
-1 0 1 2 3
0 0.9 0.1
case 7:
3
3
-1 0 1 2 3 8 10 12 15 20
0 0.5952381 0.38392857 0.02083333 0 0
Hints
Notice that N_{p}
depends only on N_{p-1}
. That means rather than recursing from the top-down, you could start at the base case and work up. While this is faster execution-wise, I'm not sure if this will necessarily produce the shortest code.
Here's a reference implementation using Python3+NumPy for the bottom-up approach:
from numpy import *
def basis(knots, m, u):
# base case
Ns = 1.*logical_and(knots[:-1] <= u, u < knots[1:])
# work backwards, starting at the base case up to p=m
for p in range(1,m+1):
Ns = (u - knots[:-1-p]) * Ns[:-1] / (knots[p:-1] - knots[:-p-1]) \
+ (knots[p+1:] - u) * Ns[1:] / (knots[p+1:]- knots[1:-p])
return Ns
Scoring
This is code golf; shortest answer in bytes wins. Standard loopholes are prohibited. You are welcome to use any built-ins desired.
It might be worth noting that the B-splines are equivalent to regular splines. – flawr – 2016-11-26T22:22:01.623
the use of
u
as both the parameter and as a knot variable was very confusing to me for a few minutes before I realized they were unrelated since I am completely knew to this stuff and have no clue what any of it actually does. Maybe consider using different variables for that? or is it standard in the literature? – Maltysen – 2016-11-26T22:30:35.440also are the knots one or zero indexed? – Maltysen – 2016-11-26T22:32:46.937
What is
k_0,k_1,...
in your examples? – flawr – 2016-11-26T22:50:17.347@flawr oops, that's a typo. Those are suppose to be the knots
u_0, u_1, ...
– helloworld922 – 2016-11-26T22:51:10.110@Maltysen I just copied the notation of Dr. Shene. It's not necessarily standard notation. I have the knots being zero-indexed, but I think the formulas still work when one-indexed (in that case, the basis functions would also be 1-indexed). – helloworld922 – 2016-11-26T22:53:26.627
@helloworld922 I suggest swapping the
– flawr – 2016-11-26T23:14:07.257non-descending
with astrictly ascending
, otherwise you could as well just remove all duplicates first to produce a meaningfull spline. (Results in the same as considering the division by 0 as 0). (Id consider this as unnecessary fluff)@flawr oh interesting, I had no idea that was a property of B-splines! I always hated the duplicate knots corner cases in my own codes. I'll modify the question accordingly. – helloworld922 – 2016-11-26T23:48:30.023
May we always take the input of
u
as a floating-point value even when it is supposed to be an integer? – R. Kap – 2016-11-27T07:07:37.417@R.Kap do you mean
m
?u
is already a real number. The actual type ofm
is not too important, so long as it can represent non-negative small integers. – helloworld922 – 2016-11-27T07:46:50.320@helloworld922 I mean even if
u=1
, can the input be taken as1.0
instead of just1
? I'm asking because in Python 2, integer division is default and to activate float division, either the numerator or denominator must be a floating-point value. – R. Kap – 2016-11-27T07:49:34.230oh, yes. That's fine. – helloworld922 – 2016-11-27T07:50:02.343