dfspvd Subroutine

private subroutine dfspvd(t, k, x, ileft, vnikx, nderiv)

Calculates value and derivs of all B-splines which do not vanish at X

Fill VNIKX(J,IDERIV), J=IDERIV, ... ,K with nonzero values of B-splines of order K+1-IDERIV , IDERIV=NDERIV, ... ,1, by repeated calls to DFSPVN

Revision history

  • 780801 DATE WRITTEN
  • 890531 Changed all specific intrinsics to generic. (WRB)
  • 890831 Modified array declarations. (WRB)
  • 890911 Removed unnecessary intrinsics. (WRB)
  • 891214 Prologue converted to Version 4.0 format. (BAB)
  • 900328 Added TYPE section. (WRB)

Arguments

Type IntentOptional Attributes Name
real(kind=wp) :: t(*)
integer :: k
real(kind=wp) :: x
integer :: ileft
real(kind=wp) :: vnikx(k,*)
integer :: nderiv

Calls

proc~~dfspvd~~CallsGraph proc~dfspvd bspline_defc_module::dfspvd proc~dfspvn bspline_defc_module::dfspvn proc~dfspvd->proc~dfspvn

Called by

proc~~dfspvd~~CalledByGraph proc~dfspvd bspline_defc_module::dfspvd proc~dfcmn bspline_defc_module::dfcmn proc~dfcmn->proc~dfspvd proc~dfc bspline_defc_module::dfc proc~dfc->proc~dfcmn

Source Code

   subroutine dfspvd (t, k, x, ileft, vnikx, nderiv)

   real(wp) :: t(*)
   integer :: k
   real(wp) :: x
   integer :: ileft
   real(wp) :: vnikx(k,*)
   integer :: nderiv

   real(wp) :: a(20,20)
   integer :: ideriv,idervm,i,j,kmd,m,jm1,ipkmd,l,jlow
   real(wp) :: fkmd,diff,v
   integer :: dfspvn_j
   real(wp), dimension(20) :: dfspvn_deltam, dfspvn_deltap

   ! set up variables for dfspvn
   dfspvn_j = 1
   dfspvn_deltam = 0.0_wp
   dfspvn_deltap = 0.0_wp

   call dfspvn(t,k+1-nderiv,1,x,ileft,vnikx(nderiv,nderiv),&
               dfspvn_j,dfspvn_deltam,dfspvn_deltap)
   if (nderiv <= 1) return

   ideriv = nderiv
   do i=2,nderiv
      idervm = ideriv-1
      do j=ideriv,k
          vnikx(j-1,idervm) = vnikx(j,ideriv)
      end do
      ideriv = idervm
      call dfspvn(t,0,2,x,ileft,vnikx(ideriv,ideriv),&
                  dfspvn_j,dfspvn_deltam,dfspvn_deltap)
   end do

   do i=1,k
      do j=1,k
         a(i,j) = 0.0_wp
      end do
      a(i,i) = 1.0_wp
   end do
   kmd = k
   do m=2,nderiv
      kmd = kmd-1
      fkmd = kmd
      i = ileft
      j = k
      do
         jm1 = j-1
         ipkmd = i + kmd
         diff = t(ipkmd) - t(i)
         if (jm1 == 0) exit
         if (diff /= 0.0_wp) then
            do l=1,j
               a(l,j) = (a(l,j) - a(l,j-1))/diff*fkmd
            end do
         end if
         j = jm1
         i = i - 1
      end do
      if (diff /= 0.0_wp) then
         a(1,1) = a(1,1)/diff*fkmd
      end if
      do i=1,k
         v = 0.0_wp
         jlow = max(i,m)
         do j=jlow,k
            v = a(i,j)*vnikx(j,m) + v
         end do
         vnikx(i,m) = v
      end do
   end do

   end subroutine dfspvd