s_bsplines_non_uniform__eval_deriv Subroutine

  SLL_PURE subroutine s_bsplines_non_uniform__eval_deriv( self, x, derivs, jmin )
    class(sll_t_bsplines_non_uniform), intent(in   ) :: self
    real(wp)                         , intent(in   ) :: x
    real(wp)                         , intent(  out) :: derivs(0:)
    integer                          , intent(  out) :: jmin

    integer  :: icell

    integer  :: j, r
    real(wp) :: temp
    real(wp) :: saved

    ! GFortran: to allocate on stack use -fstack-arrays
    real(wp) :: left (1:self%degree)
    real(wp) :: right(1:self%degree)

    ! Check on inputs
    SLL_ASSERT( x > self%xmin - 1.0d-14 )
    SLL_ASSERT( x < self%xmax + 1.0d-14 )
    SLL_ASSERT( size(derivs) == 1+self%degree )

    ! 1. Compute cell index 'icell' and x_offset
!    icell = f_bsplines_non_uniform__find_cell( self, x )
    icell = self % find_cell( x )
    SLL_ASSERT( icell >= 1               )
    SLL_ASSERT( icell <= self%ncells     )
    SLL_ASSERT( self%knots(icell) <= x   )
    SLL_ASSERT( x <= self%knots(icell+1) )

    ! 2. Compute index range of B-splines with support over cell 'icell'
    jmin = icell - self%offset

    ! 3. Compute derivatives of aforementioned B-splines
    associate( degree => self%degree, degree_real => real(self%degree,wp) )

      ! Compute nonzero basis functions and knot differences
      ! for splines up to degree deg-1 which are needed to compute derivative
      ! First part of Algorithm  A3.2 of NURBS book
      derivs(0) = 1.0_wp
      do j = 1, degree-1
         left (j) = x - self%knots(icell+1-j)
         right(j) = self%knots(icell+j) - x
         saved    = 0.0_wp
         do r = 0, j-1
            ! compute and save bspline values
            temp      = derivs(r)/(right(r+1) + left(j-r))
            derivs(r) = saved + right(r+1) * temp
            saved     = left(j-r) * temp
         end do
         derivs(j) = saved
      end do

      ! Compute derivatives at x using values stored in bsdx and formula
      ! formula for spline derivative based on difference of splines of
      ! degree deg-1
      ! -------
      ! j = 0
      saved = degree_real * derivs(0) / (self%knots(icell+1)-self%knots(icell+1-degree))
      derivs(0) = -saved
      do j = 1, degree-1
         temp      = saved
         saved     = degree_real * derivs(j) / (self%knots(icell+j+1)-self%knots(icell+j+1-degree))
         derivs(j) = temp - saved
      end do
      ! j = degree
      derivs(degree) =  saved

    end associate

  end subroutine s_bsplines_non_uniform__eval_deriv