code: def find_ctrlpts_surface(t_u, t_v, surf, **kwargs):
      
      # Get keyword arguments
      span_func = kwargs.get('find_span_func', helpers.find_span_linear)
  
      # Find spans
      span_u = span_func(surf.degree_u, surf.knotvector_u, surf.ctrlpts_size_u, t_u)
      span_v = span_func(surf.degree_v, surf.knotvector_v, surf.ctrlpts_size_v, t_v)
  
      # Constant indices
      idx_u = span_u - surf.degree_u
      idx_v = span_v - surf.degree_v
  
      # Find control points involved in evaluation of the surface point at the input parameter pair (u, v)
      surf_ctrlpts = [[] for _ in range(surf.degree_u + 1)]
      for k in range(surf.degree_u + 1):
          temp = [() for _ in range(surf.degree_v + 1)]
          for l in range(surf.degree_v + 1):
              temp[l] = surf.ctrlpts2d[idx_u + k][idx_v + l]
          surf_ctrlpts[k] = temp
  
      # Return 2-dimensional control points array
      return surf_ctrlpts
docs: Finds the control points involved in the evaluation of the surface point defined by the input parameter pair.
  
      This function uses a modified version of the algorithm *A3.5 SurfacePoint* from The NURBS Book by Piegl & Tiller.
  
      :param t_u: parameter on the u-direction
      :type t_u: float
      :param t_v: parameter on the v-direction
      :type t_v: float
      :param surf: input surface
      :type surf: abstract.Surface
      :return: 2-dimensional control points array
      :rtype: list