dbtpcf computes b-spline interpolation coefficients for nf sets of data stored in the columns of the array fcn. the b-spline coefficients are stored in the rows of bcoef however. each interpolation is based on the n abcissa stored in the array x, and the n+k knots stored in the array t. the order of each interpolation is k.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real(kind=wp), | intent(in), | dimension(:) | :: | x | ||
integer(kind=ip), | intent(in) | :: | n |
dimension of |
||
real(kind=wp), | intent(in), | dimension(ldf,nf) | :: | fcn | ||
integer(kind=ip), | intent(in) | :: | ldf | |||
integer(kind=ip), | intent(in) | :: | nf | |||
real(kind=wp), | intent(in), | dimension(:) | :: | t | ||
integer(kind=ip), | intent(in) | :: | k | |||
real(kind=wp), | intent(out), | dimension(nf,n) | :: | bcoef | ||
real(kind=wp), | intent(out), | dimension(*) | :: | work |
work array of size >= |
|
integer(kind=ip), | intent(out) | :: | iflag |
status flag:
|
pure subroutine dbtpcf(x,n,fcn,ldf,nf,t,k,bcoef,work,iflag) integer(ip),intent(in) :: n !! dimension of `x` integer(ip),intent(in) :: nf integer(ip),intent(in) :: ldf integer(ip),intent(in) :: k real(wp),dimension(:),intent(in) :: x real(wp),dimension(ldf,nf),intent(in) :: fcn real(wp),dimension(:),intent(in) :: t real(wp),dimension(nf,n),intent(out) :: bcoef real(wp),dimension(*),intent(out) :: work !! work array of size >= `2*k*(n+1)` integer(ip),intent(out) :: iflag !! status flag: !! !! * 0: no errors !! * 301: n should be >0 integer(ip) :: i, j, m1, m2, iq, iw ! check for null input if (nf > 0_ip) then ! partition work array m1 = k - 1_ip m2 = m1 + k iq = 1_ip + n iw = iq + m2*n+1_ip ! compute b-spline coefficients ! first data set call dbintk(x,fcn,t,n,k,work,work(iq),work(iw),iflag) if (iflag == 0_ip) then do i=1_ip,n bcoef(1_ip,i) = work(i) end do ! all remaining data sets by back-substitution if (nf == 1_ip) return do j=2_ip,nf do i=1_ip,n work(i) = fcn(i,j) end do call dbnslv(work(iq),m2,n,m1,m1,work) do i=1_ip,n bcoef(j,i) = work(i) end do end do end if else !write(error_unit,'(A)') 'dbtpcf - n should be >0' iflag = 301_ip end if end subroutine dbtpcf