s_spline_interpolator_1d__compute_interpolant Subroutine

  subroutine s_spline_interpolator_1d__compute_interpolant( self, &
      bcoef, gtau, derivs_xmin, derivs_xmax )
    class(sll_t_spline_interpolator_1d), intent(in   ) :: self
  !  type (sll_t_spline_1d)             , intent(inout) :: spline
    real(wp)                           , intent(  out) :: &  ! not using spline to separate data and basis
    bcoef(1-MERGE(1+self%bspl%degree/2,0,self%bspl%periodic):self%bspl%nbasis+MERGE(1+self%bspl%degree/2,0,self%bspl%periodic))
    real(wp)                           , intent(in   ) :: gtau(:)
    real(wp),                  optional, intent(in   ) :: derivs_xmin(:)
    real(wp),                  optional, intent(in   ) :: derivs_xmax(:)

    character(len=*), parameter :: &
      this_sub_name = "sll_t_spline_interpolator_1d % compute_interpolant"

    integer :: i

    SLL_ASSERT( size(gtau) == self%bspl%nbasis - self%nbc_xmin - self%nbc_xmax )
!    SLL_ASSERT( spline % belongs_to_space( self % bspl ) )

    associate ( nbasis   =>   self % bspl % nbasis, &
                degree   =>   self % bspl % degree, &
                nbc_xmin =>   self % nbc_xmin     , &
                nbc_xmax =>   self % nbc_xmax     )!, &
!                bcoef    => spline % bcoef )

      ! Special case: linear spline
      if (degree == 1) then
        bcoef(1:nbasis) = gtau(1:nbasis)
        ! Periodic only: "wrap around" coefficients onto extended array
        if (self%bspl%periodic) then
          bcoef(0)        = bcoef(nbasis)
          bcoef(nbasis+1) = bcoef(1)
        end if
        return
      end if

      ! Hermite boundary conditions at xmin, if any
      ! NOTE: For consistency with the linear system, the i-th derivative
      !       provided by the user must be multiplied by dx^i
      if ( self%bc_xmin == sll_p_hermite ) then
        if ( present( derivs_xmin ) ) then
          bcoef(1:nbc_xmin) = &
                     [(derivs_xmin(i)*self%dx**(i+self%odd-1), i=nbc_xmin,1,-1)]
        else
          SLL_ERROR(this_sub_name,"Hermite BC at xmin requires derivatives")
        end if
      end if

      ! Interpolation points
      bcoef(nbc_xmin+1:nbasis-nbc_xmax) = gtau(:)

      ! Hermite boundary conditions at xmax, if any
      ! NOTE: For consistency with the linear system, the i-th derivative
      !       provided by the user must be multiplied by dx^i
      if ( self%bc_xmax == sll_p_hermite ) then
        if ( present( derivs_xmax ) ) then
          bcoef(nbasis-nbc_xmax+1:nbasis) = &
                        [(derivs_xmax(i)*self%dx**(i+self%odd-1), i=1,nbc_xmax)]
        else
          SLL_ERROR(this_sub_name,"Hermite BC at xmax requires derivatives")
        end if
      end if

      ! Solve linear system and compute coefficients
      call self % matrix % solve_inplace(1, bcoef(1:nbasis) )

      ! Periodic only: "wrap around" coefficients onto extended array
      if (self%bc_xmin == sll_p_periodic) then
        associate( g => 1+degree/2 )
          bcoef(1-g:0)             = bcoef(nbasis-g+1:nbasis)
          bcoef(nbasis+1:nbasis+g) = bcoef(1:g)
        end associate
      end if

    end associate ! nbasis, degree, bcoef

  end subroutine s_spline_interpolator_1d__compute_interpolant