code: def _make_interp_spline_full_matrix(x, y, k, t, bc_type):
      """Construct the interpolating spline spl(x) = y with *full* linalg.
  
      Only useful for testing, do not call directly!
      This version is O(N**2) in memory and O(N**3) in flop count.
      """
      # convert string aliases for the boundary conditions
      if bc_type is None or bc_type == "not-a-knot":
          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
  
      # 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
      # 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
      assert nt - n == nleft + nright
  
      # 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.
  
      # 1. Compute intervals for each value
      intervals = cupy.empty_like(x, dtype=cupy.int64)
      interval_kernel = _get_module_func(
          INTERVAL_MODULE, "find_interval"
      )
      interval_kernel(
          ((x.shape[0] + 128 - 1) // 128,),
          (128,),
          (t, x, intervals, k, nt, False, x.shape[0]),
      )
  
      # 2. Compute non-zero b-spline basis elements for each value in `x`
      # The way de_Boor_D kernel is written, it wants `c` and `out` arrays
      # which we do not use (but need to provide to the kernel), and the
      # `temp` array contains non-zero b-spline basis elements, which we do want.
      dummy_c = cupy.empty((nt, 1), dtype=float)
      out = cupy.empty(
          (len(x), prod(dummy_c.shape[1:])), dtype=dummy_c.dtype
      )
  
      num_c = prod(dummy_c.shape[1:])
      temp = cupy.empty(x.shape[0] * (2 * k + 1))
      d_boor_kernel = _get_module_func(D_BOOR_MODULE, "d_boor", dummy_c)
      d_boor_kernel(
          ((x.shape[0] + 128 - 1) // 128,),
          (128,),
          (
              t,
              dummy_c,
              k,
              0,
              x,
              intervals,
              out,
              temp,
              num_c,
              0,
              x.shape[0],
          ),
      )
  
      # 3. Construct the colocation matrix.
      # For each value in `x`, the `temp` array contains 2k+1 entries : first
      # k+1 elements are b-splines, followed by k entries used for work storage
      # which we ignore.
      # XXX: full matrices! Can / should use banded linear algebra instead.
      A = cupy.zeros((nt, nt), dtype=float)
      offset = nleft
      for j in range(len(x)):
          left = intervals[j]
          A[j + offset, left - k : left + 1] = temp[
              j * (2 * k + 1) : j * (2 * k + 1) + k + 1
          ]
  
      # 4. Handle boundary conditions: The colocation matrix is augmented with
      # additional rows, one row per derivative at the left and right edges.
      # We need a python loop for each derivative because in general they can be
      # of any order, `m`.
      # The left-side boundary conditions go to the first rows of the matrix
      # and the right-side boundary conditions go to the last rows.
      intervals_bc = cupy.empty(1, dtype=cupy.int64)
      if nleft > 0:
          intervals_bc[0] = intervals[0]
          x0 = cupy.array([x[0]], dtype=x.dtype)
          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,
                      temp,
                      num_c,
                      0,
                      1,
                  ),
              )
              left = intervals_bc[0]
              A[j, left - k : left + 1] = temp[: k + 1]
  
      # repeat for the b.c. at the right edge.
      if nright > 0:
          intervals_bc[0] = intervals[-1]
          x0 = cupy.array([x[-1]], dtype=x.dtype)
          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,
                      temp,
                      num_c,
                      0,
                      1,
                  ),
              )
              left = intervals_bc[0]
              row = nleft + len(x) + j
              A[row, left - k : left + 1] = temp[: k + 1]
  
      # 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 for the coefficients.
      from cupy.linalg import solve
  
      coef = solve(A, 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 10 for calling a global function with example usage: _make_interp_spline_full_matrix(x, y, k, t, bc_type) and returns: BSpline