cupyx.scipy.interpolate.BSpline#

class cupyx.scipy.interpolate.BSpline(t, c, k, extrapolate=True, axis=0)[source]#

Univariate spline in the B-spline basis.

\[S(x) = \sum_{j=0}^{n-1} c_j B_{j, k; t}(x)\]

where \(B_{j, k; t}\) are B-spline basis functions of degree k and knots t.

Parameters:
  • t (ndarray, shape (n+k+1,)) – knots

  • c (ndarray, shape (>=n, ...)) – spline coefficients

  • k (int) – B-spline degree

  • extrapolate (bool or 'periodic', optional) – whether to extrapolate beyond the base interval, t[k] .. t[n], or to return nans. If True, extrapolates the first and last polynomial pieces of b-spline functions active on the base interval. If ‘periodic’, periodic extrapolation is used. Default is True.

  • axis (int, optional) – Interpolation axis. Default is zero.

Variables:
  • t (ndarray) – knot vector

  • c (ndarray) – spline coefficients

  • k (int) – spline degree

  • extrapolate (bool) – If True, extrapolates the first and last polynomial pieces of b-spline functions active on the base interval.

  • axis (int) – Interpolation axis.

  • tck (tuple) – A read-only equivalent of (self.t, self.c, self.k)

Notes

B-spline basis elements are defined via

\[ \begin{align}\begin{aligned}B_{i, 0}(x) = 1, \textrm{if $t_i \le x < t_{i+1}$, otherwise $0$,}\\B_{i, k}(x) = \frac{x - t_i}{t_{i+k} - t_i} B_{i, k-1}(x) + \frac{t_{i+k+1} - x}{t_{i+k+1} - t_{i+1}} B_{i+1, k-1}(x)\end{aligned}\end{align} \]

Implementation details

  • At least k+1 coefficients are required for a spline of degree k, so that n >= k+1. Additional coefficients, c[j] with j > n, are ignored.

  • B-spline basis elements of degree k form a partition of unity on the base interval, t[k] <= x <= t[n].

  • Based on [1] and [2]

References

Methods

__call__(x, nu=0, extrapolate=None)[source]#

Evaluate a spline function.

Parameters:
  • x (array_like) – points to evaluate the spline at.

  • nu (int, optional) – derivative to evaluate (default is 0).

  • extrapolate (bool or 'periodic', optional) – whether to extrapolate based on the first and last intervals or return nans. If ‘periodic’, periodic extrapolation is used. Default is self.extrapolate.

Returns:

y – Shape is determined by replacing the interpolation axis in the coefficient array with the shape of x.

Return type:

array_like

antiderivative(nu=1)[source]#

Return a B-spline representing the antiderivative.

Parameters:

nu (int, optional) – Antiderivative order. Default is 1.

Returns:

b – A new instance representing the antiderivative.

Return type:

BSpline object

Notes

If antiderivative is computed and self.extrapolate='periodic', it will be set to False for the returned instance. This is done because the antiderivative is no longer periodic and its correct evaluation outside of the initially given x interval is difficult.

See also

splder, splantider

classmethod basis_element(t, extrapolate=True)[source]#

Return a B-spline basis element B(x | t[0], ..., t[k+1]).

Parameters:
  • t (ndarray, shape (k+2,)) – internal knots

  • extrapolate (bool or 'periodic', optional) – whether to extrapolate beyond the base interval, t[0] .. t[k+1], or to return nans. If ‘periodic’, periodic extrapolation is used. Default is True.

Returns:

basis_element – A callable representing a B-spline basis element for the knot vector t.

Return type:

callable

Notes

The degree of the B-spline, k, is inferred from the length of t as len(t)-2. The knot vector is constructed by appending and prepending k+1 elements to internal knots t.

classmethod construct_fast(t, c, k, extrapolate=True, axis=0)[source]#

Construct a spline without making checks. Accepts same parameters as the regular constructor. Input arrays t and c must of correct shape and dtype.

derivative(nu=1)[source]#

Return a B-spline representing the derivative.

Parameters:

nu (int, optional) – Derivative order. Default is 1.

Returns:

b – A new instance representing the derivative.

Return type:

BSpline object

See also

splder, splantider

classmethod design_matrix(x, t, k, extrapolate=False)[source]#

Returns a design matrix as a CSR format sparse array.

Parameters:
  • x (array_like, shape (n,)) – Points to evaluate the spline at.

  • t (array_like, shape (nt,)) – Sorted 1D array of knots.

  • k (int) – B-spline degree.

  • extrapolate (bool or 'periodic', optional) – Whether to extrapolate based on the first and last intervals or raise an error. If ‘periodic’, periodic extrapolation is used. Default is False.

Returns:

design_matrix – Sparse matrix in CSR format where each row contains all the basis elements of the input row (first row = basis elements of x[0], …, last row = basis elements x[-1]).

Return type:

csr_matrix object

Notes

In each row of the design matrix all the basis elements are evaluated at the certain point (first row - x[0], …, last row - x[-1]). nt is a length of the vector of knots: as far as there are nt - k - 1 basis elements, nt should be not less than 2 * k + 2 to have at least k + 1 basis element.

Out of bounds x raises a ValueError.

Note

This method returns a csr_matrix instance as CuPy still does not have csr_array.

integrate(a, b, extrapolate=None)[source]#

Compute a definite integral of the spline.

Parameters:
  • a (float) – Lower limit of integration.

  • b (float) – Upper limit of integration.

  • extrapolate (bool or 'periodic', optional) – whether to extrapolate beyond the base interval, t[k] .. t[-k-1], or take the spline to be zero outside of the base interval. If ‘periodic’, periodic extrapolation is used. If None (default), use self.extrapolate.

Returns:

I – Definite integral of the spline over the interval [a, b].

Return type:

array_like

__eq__(value, /)#

Return self==value.

__ne__(value, /)#

Return self!=value.

__lt__(value, /)#

Return self<value.

__le__(value, /)#

Return self<=value.

__gt__(value, /)#

Return self>value.

__ge__(value, /)#

Return self>=value.

Attributes

tck#

Equivalent to (self.t, self.c, self.k) (read-only).