dpchkt Subroutine

private subroutine dpchkt(n, x, knotyp, t)

Set a knot sequence for the B-spline representation of a PCH function with breakpoints X. All knots will be at least double. Endknots are set as:

  1. quadruple knots at endpoints if KNOTYP=0;
  2. extrapolate the length of end interval if KNOTYP=1;
  3. periodic if KNOTYP=2.

Restrictions/assumptions:

  1. N>=2 . (not checked)
  2. X(i)<X(i+1), i=1,...,N . (not checked)
  3. 0<=KNOTYP<=2 . (Acts like KNOTYP=0 for any other value.)

Internal Notes

Since this is subsidiary to dpchbs, which validates its input before calling, it is unnecessary for such validation to be done here.

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: n
real(kind=wp), intent(in) :: x(*)
integer, intent(in) :: knotyp
real(kind=wp), intent(out) :: t(*)

Called by

proc~~dpchkt~~CalledByGraph proc~dpchkt pchip_module::dpchkt proc~dpchbs pchip_module::dpchbs proc~dpchbs->proc~dpchkt

Source Code

    subroutine dpchkt (n, x, knotyp, t)

    integer,intent(in)   :: n
    integer,intent(in)   :: knotyp
    real(wp),intent(in)  :: x(*)
    real(wp),intent(out) :: t(*)

    integer :: j, k, ndim
    real(wp) :: hbeg, hend

    ! initialize.
    ndim = 2*n

    ! set interior knots.
    j = 1
    do k = 1, n
        j = j + 2
        t(j) = x(k)
        t(j+1) = t(j)
    end do

    ! assertion:  at this point t(3),...,t(ndim+2) have been set and
    !             j=ndim+1.
    !
    !  set end knots according to knotyp.

    hbeg = x(2) - x(1)
    hend = x(n) - x(n-1)
    select case (knotyp)
    case (1)
        ! extrapolate.
        t(2) = x(1) - hbeg
        t(ndim+3) = x(n) + hend
    case (2)
        ! periodic.
        t(2) = x(1) - hend
        t(ndim+3) = x(n) + hbeg
    case default
        ! quadruple end knots.
        t(2) = x(1)
        t(ndim+3) = x(n)
    end select
    t(1) = t(2)
    t(ndim+4) = t(ndim+3)

    end subroutine dpchkt