code: def make_interp_spline(
      x, y, k=3, t=None, bc_type=None, axis=0, check_finite=True
  ):
      """Compute the (coefficients of) interpolating B-spline.
  
      Parameters
      ----------
      x : array_like, shape (n,)
          Abscissas.
      y : array_like, shape (n, ...)
          Ordinates.
      k : int, optional
          B-spline degree. Default is cubic, ``k = 3``.
      t : array_like, shape (nt + k + 1,), optional.
          Knots.
          The number of knots needs to agree with the number of data points and
          the number of derivatives at the edges. Specifically, ``nt - n`` must
          equal ``len(deriv_l) + len(deriv_r)``.
      bc_type : 2-tuple or None
          Boundary conditions.
          Default is None, which means choosing the boundary conditions
          automatically. Otherwise, it must be a length-two tuple where the first
          element (``deriv_l``) sets the boundary conditions at ``x[0]`` and
          the second element (``deriv_r``) sets the boundary conditions at
          ``x[-1]``. Each of these must be an iterable of pairs
          ``(order, value)`` which gives the values of derivatives of specified
          orders at the given edge of the interpolation interval.
          Alternatively, the following string aliases are recognized:
  
          * ``"clamped"``: The first derivatives at the ends are zero. This is
             equivalent to ``bc_type=([(1, 0.0)], [(1, 0.0)])``.
          * ``"natural"``: The second derivatives at ends are zero. This is
            equivalent to ``bc_type=([(2, 0.0)], [(2, 0.0)])``.
          * ``"not-a-knot"`` (default): The first and second segments are the
            same polynomial. This is equivalent to having ``bc_type=None``.
          * ``"periodic"``: The values and the first ``k-1`` derivatives at the
            ends are equivalent.
  
      axis : int, optional
          Interpolation axis. Default is 0.
      check_finite : bool, optional
          Whether to check that the input arrays contain only finite numbers.
          Disabling may give a performance gain, but may result in problems
          (crashes, non-termination) if the inputs do contain infinities or NaNs.
          Default is True.
  
      Returns
      -------
      b : a BSpline object of the degree ``k`` and with knots ``t``.
  
      """
      # convert string aliases for the boundary conditions
      if (
          bc_type is None
          or bc_type == "not-a-knot"
          or bc_type == "periodic"
      ):
          deriv_l, deriv_r = None, None
      elif isinstance(bc_type, str):
          deriv_l, deriv_r = bc_type, bc_type
      else:
          try:
              deriv_l, deriv_r = bc_type
          except TypeError as e:
              raise ValueError(
                  "Unknown boundary condition: %s" % bc_type
              ) from e
  
      y = cupy.asarray(y)
  
      axis = normalize_axis_index(axis, y.ndim)
  
      x = _as_float_array(x, check_finite)
      y = _as_float_array(y, check_finite)
  
      y = cupy.moveaxis(
          y, axis, 0
      )  # now internally interp axis is zero
  
      # sanity check the input
      if bc_type == "periodic" and not cupy.allclose(
          y[0], y[-1], atol=1e-15
      ):
          raise ValueError(
              "First and last points does not match while "
              "periodic case expected"
          )
      if x.size != y.shape[0]:
          raise ValueError(
              "Shapes of x {} and y {} are incompatible".format(
                  x.shape, y.shape
              )
          )
      if (x[1:] == x[:-1]).any():
          raise ValueError("Expect x to not have duplicates")
      if x.ndim != 1 or (x[1:] < x[:-1]).any():
          raise ValueError(
              "Expect x to be a 1D strictly increasing sequence."
          )
  
      # special-case k=0 right away
      if k == 0:
          if any(_ is not None for _ in (t, deriv_l, deriv_r)):
              raise ValueError(
                  "Too much info for k=0: t and bc_type can only "
                  "be None."
              )
          t = cupy.r_[x, x[-1]]
          c = cupy.asarray(y)
          c = cupy.ascontiguousarray(c, dtype=_get_dtype(c.dtype))
          return BSpline.construct_fast(t, c, k, axis=axis)
  
      # special-case k=1 (e.g., Lyche and Morken, Eq.(2.16))
      if k == 1 and t is None:
          if not (deriv_l is None and deriv_r is None):
              raise ValueError(
                  "Too much info for k=1: bc_type can only be None."
              )
          t = cupy.r_[x[0], x, x[-1]]
          c = cupy.asarray(y)
          c = cupy.ascontiguousarray(c, dtype=_get_dtype(c.dtype))
          return BSpline.construct_fast(t, c, k, axis=axis)
  
      k = operator.index(k)
  
      if bc_type == "periodic" and t is not None:
          raise NotImplementedError(
              "For periodic case t is constructed "
              "automatically and can not be passed "
              "manually"
          )
  
      # come up with a sensible knot vector, if needed
      if t is None:
          if deriv_l is None and deriv_r is None:
              if bc_type == "periodic":
                  t = _periodic_knots(x, k)
              elif k == 2:
                  # OK, it's a bit ad hoc: Greville sites + omit
                  # 2nd and 2nd-to-last points, a la not-a-knot
                  t = (x[1:] + x[:-1]) / 2.0
                  t = cupy.r_[
                      (x[0],) * (k + 1), t[1:-1], (x[-1],) * (k + 1)
                  ]
              else:
                  t = _not_a_knot(x, k)
          else:
              t = _augknt(x, k)
  
      t = _as_float_array(t, check_finite)
  
      if k < 0:
          raise ValueError("Expect non-negative k.")
      if t.ndim != 1 or (t[1:] < t[:-1]).any():
          raise ValueError("Expect t to be a 1-D sorted array_like.")
      if t.size < x.size + k + 1:
          raise ValueError(
              "Got %d knots, need at least %d."
              % (t.size, x.size + k + 1)
          )
      if (x[0] < t[k]) or (x[-1] > t[-k]):
          raise ValueError("Out of bounds w/ x = %s." % x)
  
      if bc_type == "periodic":
          return _make_periodic_spline(x, y, t, k, axis)
  
      # Here : deriv_l, r = [(nu, value), ...]
      deriv_l = _convert_string_aliases(deriv_l, y.shape[1:])
      deriv_l_ords, deriv_l_vals = _process_deriv_spec(deriv_l)
      nleft = deriv_l_ords.shape[0]
  
      deriv_r = _convert_string_aliases(deriv_r, y.shape[1:])
      deriv_r_ords, deriv_r_vals = _process_deriv_spec(deriv_r)
      nright = deriv_r_ords.shape[0]
  
      # have `n` conditions for `nt` coefficients; need nt-n derivatives
      n = x.size
      nt = t.size - k - 1
  
      if nt - n != nleft + nright:
          raise ValueError(
              "The number of derivatives at boundaries does not "
              "match: expected %s, got %s + %s"
              % (nt - n, nleft, nright)
          )
  
      # bail out if the `y` array is zero-sized
      if y.size == 0:
          c = cupy.zeros((nt,) + y.shape[1:], dtype=float)
          return BSpline.construct_fast(t, c, k, axis=axis)
  
      # Consruct the colocation matrix of b-splines + boundary conditions.
      # The coefficients of the interpolating B-spline function are the solution
      # of the linear system `A @ c = rhs` where `A` is the colocation matrix
      # (i.e., each row of A corresponds to a data point in the `x` array and
      #  contains b-splines which are non-zero at this value of x)
      # Each boundary condition is a fixed value of a certain derivative
      # at the edge, so each derivative adds a row to `A`.
      # The `rhs` is the array of data values, `y`, plus derivatives from
      # boundary conditions, if any.
      # The colocation matrix is banded (has at most k+1 diagonals). Since LAPACK
      # linear algebra (?gbsv) is not available, we store it as a CSR array
  
      # 1. Construct the colocation matrix itself.
      matr = BSpline.design_matrix(x, t, k)
  
      # 2. Boundary conditions: need to augment the design matrix with additional
      # rows, one row per derivative at the left and right edges.
      # The left-side boundary conditions go to the first rows of the matrix
      # and the right-side boundary conditions go to the last rows.
      # Will need a python loop for each derivative because in general they
      # can be of any order, `m`.
      # To compute the derivatives, will invoke the de Boor D kernel.
      if nleft > 0 or nright > 0:
          # Prepare the I/O arrays for the kernels. We only need the non-zero
          # b-splines at x[0] and x[-1], but the kernel wants more arrays which
          # we allocate and ignore (mode != 1)
          temp = cupy.zeros((nt,), dtype=float)
          num_c = 1
          dummy_c = cupy.empty((nt, num_c), dtype=float)
          out = cupy.empty((1, 1), dtype=dummy_c.dtype)
  
          d_boor_kernel = _get_module_func(
              D_BOOR_MODULE, "d_boor", dummy_c
          )
  
          # find the intervals for x[0] and x[-1]
          intervals_bc = cupy.empty(2, dtype=cupy.int64)
          interval_kernel = _get_module_func(
              INTERVAL_MODULE, "find_interval"
          )
          interval_kernel(
              (1,),
              (2,),
              (t, cupy.r_[x[0], x[-1]], intervals_bc, k, nt, False, 2),
          )
  
      # 3. B.C.s at x[0]
      if nleft > 0:
          x0 = cupy.array([x[0]], dtype=x.dtype)
          rows = cupy.zeros((nleft, nt), dtype=float)
  
          for j, m in enumerate(deriv_l_ords):
              # place the derivatives of the order m at x[0] into `temp`
              d_boor_kernel(
                  (1,),
                  (1,),
                  (
                      t,
                      dummy_c,
                      k,
                      int(m),
                      x0,
                      intervals_bc,
                      out,  # ignore (mode !=1),
                      temp,  # non-zero b-splines
                      num_c,  # the 2nd dimension of `dummy_c`. Ignore.
                      0,  # mode != 1 => do not touch dummy_c array
                      1,
                  ),
              )  # the length of the `x0` array
              left = intervals_bc[0]
              rows[j, left - k : left + 1] = temp[: k + 1]
  
          matr = sparse.vstack(
              [sparse.csr_matrix(rows), matr]  # A[:nleft, :]
          )
  
      # 4. Repeat for B.C.s at x[-1]
      if nright > 0:
          intervals_bc[0] = intervals_bc[
              -1
          ]  # use the intervals for x[-1]
          x0 = cupy.array([x[-1]], dtype=x.dtype)
          rows = cupy.zeros((nright, nt), dtype=float)
  
          for j, m in enumerate(deriv_r_ords):
              # place the derivatives of the order m at x[0] into `temp`
              d_boor_kernel(
                  (1,),
                  (1,),
                  (
                      t,
                      dummy_c,
                      k,
                      int(m),
                      x0,
                      intervals_bc,
                      out,  # ignore (mode !=1),
                      temp,  # non-zero b-splines
                      num_c,  # the 2nd dimension of `dummy_c`. Ignore.
                      0,  # mode != 1 => do not touch dummy_c array
                      1,
                  ),
              )  # the length of the `x0` array
              left = intervals_bc[0]
              rows[j, left - k : left + 1] = temp[: k + 1]
  
          matr = sparse.vstack(
              [matr, sparse.csr_matrix(rows)]
          )  # A[nleft+len(x):, :]
  
      # 5. Prepare the RHS: `y` values to interpolate (+ derivatives, if any)
      extradim = prod(y.shape[1:])
      rhs = cupy.empty((nt, extradim), dtype=y.dtype)
      if nleft > 0:
          rhs[:nleft] = deriv_l_vals.reshape(-1, extradim)
      rhs[nleft : nt - nright] = y.reshape(-1, extradim)
      if nright > 0:
          rhs[nt - nright :] = deriv_r_vals.reshape(-1, extradim)
  
      # 6. Finally, solve the linear system for the coefficients.
      if cupy.issubdtype(rhs.dtype, cupy.complexfloating):
          # avoid upcasting the l.h.s. to complex (that doubles the memory)
          coef = (
              spsolve(matr, rhs.real) + spsolve(matr, rhs.imag) * 1.0j
          )
      else:
          coef = spsolve(matr, rhs)
      coef = cupy.ascontiguousarray(coef.reshape((nt,) + y.shape[1:]))
      return BSpline(t, coef, k)
  
file_path: cupy/cupyx/scipy/interpolate/_bspline2.py
imports: import operator
  from numpy.core.multiarray import normalize_axis_index
  import cupy
  from cupyx.scipy import sparse
  from cupyx.scipy.sparse.linalg import spsolve
  from cupyx.scipy.interpolate._bspline import _get_module_func, INTERVAL_MODULE, D_BOOR_MODULE, BSpline
label_desc: Use image node_id 9 for calling a global function with example usage: make_interp_spline(x, y, k, t, bc_type, axis, check_finite) and returns: BSpline, BSpline, BSpline, _make_periodic_spline, BSpline