◆ Pchic()

Sub Pchic	(	Ic() As	Long,
		Vc() As	Double,
		Switch As	Double,
		N As	Long,
		X() As	Double,
		F() As	Double,
		D() As	Double,
		Info As	Long,
		Optional Idxfd As	Long = `0`,
		Optional Incfd As	Long = `1`
	)

Piecewise cubic Hermite interpolation

Purpose: This routine sets derivatives needed to determine a monotone piecewise cubic Hermite interpolant. User control is available over boundary conditions and/or treatment of points where monotonicity switches direction.

The resulting piecewise cubic Hermite function may be evaluated by Pchfe or Pchfd.

Parameters

[in]	Ic()	Array Ic(LIc - 1) (LIc >= 2) Specification of desired boundary conditions. Ic(0) = Desired condition at beginning of data. = 0: The default boundary condition (the same as used by Pchim). > 0: The boundary derivative is not to be adjusted for monotonicity. < 0: The derivative is to be adjusted for monotonicity. Allowable values for the magnitude of Ic(0) are: = 1: if first derivative at X(0) is given in Vc(0) = 2: if second derivative at X(0) is given in Vc(0) = 3: to use the 3-point difference formula for D(0) (reverts to the default boundary conditions if N < 3) = 4: to use the 4-point difference formula for D(0) (reverts to the default boundary conditions if N < 4) = 5: to set D(0) so that the second derivative is continuous at X(1) (reverts to the default boundary conditions if N < 4) Ic(1) = Desired condition at end of data. It may take on the same values as Ic(0), but applied to derivative at X(N-1). In the case Ic(1) = 1 or 2, the value is given in Vc(1)
[in]	Vc()	Array Vc(LVc - 1) (LVc >= 2) Specification of desired boundary values, as indicated above. Vc(0) need to be set only if Ic(0) = 1 or 2. Vc(1) need to be set only if Ic(1) = 1 or 2.
[in]	Switch	Indicates desired treatment of points where direction of monotonicity switches. = 0: Interpolant is required to be monotonic in each interval, regardless of monotonicity of data. <> 0: Use a formula based on the 3-point difference formula in the vicinity of switch points. > 0: The interpolant on each interval containing an extremum is controlled to not deviate from the data by more than Switch * Dfloc, where Dfloc is the maximum of the change of F() on this interval and its two immediate neighbors. < 0: No such control is to be imposed.
[in]	N	Number of data points. (N >= 2)
[in]	X()	Array X(LX - 1) (LX >= N) Independent variable values. The elements of X() must be strictly increasing.
[in]	F()	Array F(LF - 1) (LF >= Incfd*(N - 1) + 1) Dependent variable values to be interpolated.
[out]	D()	Array D(LD - 1) (LD >= Incfd*(N - 1) + 1) Derivative values at the data points.
[out]	Info	= 0: Successful exit. = -1: The argument Ic() had an illegal value. (Abs(Ic(0)) > 5 and/or Abs(Ic(1)) > 5) = -2: The argument Vc() is invalid. = -4: The argument N had an illegal value. (N < 2) = -5: The argument X() had an illegal value. (X() is not strictly increasing) = -6: The argument F() is invalid. = -7: The argument D() is invalid. = -11: The argument Incfd had an illegal value. (Incfd < 1) = 1 Ic(0) < 0 and D(0) had to be adjusted for monotonicity. = 2 Ic(1) < 0 and D((N-1)*Incfd) had to be adjusted for monotonicity. = 3 Both of the above are true.
[in]	Idxfd	(Optional) Index of the first data in F() and D(). (default = 0)
[in]	Incfd	(Optional) Increment between successive values in F() and D(). (Incfd >= 1) (default = 1)

Reference: SLATEC (PCHIP)

Example Program

Compute ln(0.115) and ln(0.125) by interpolating the following natural logarithm table with cubic spline.
In this example, Pchfd is used to evaluate the interpolation function instead of Pchfe, so that the derivatives ln'(0.115) and ln((0.125) are also computed.

  x       ln(x)
------ ---------
 0.10   -2.3026
 0.11   -2.2073
 0.12   -2.1203
 0.13   -2.0402
------ ---------

Sub Ex_Pchic()
    Const N = 4, Ne = 2
    Dim X(N - 1) As Double, Y(N - 1) As Double, D(N - 1) As Double
    Dim Ic(1) As Long, Vc(1) As Double
    Dim Xe(Ne - 1) As Double, Ye(Ne - 1) As Double, Yep(Ne - 1) As Double
    Dim Info As Long, I As Long
    '-- Data
    X(0) = 0.1: Y(0) = -2.3026
    X(1) = 0.11: Y(1) = -2.2073
    X(2) = 0.12: Y(2) = -2.1203
    X(3) = 0.13: Y(3) = -2.0402
    '-- Cubic Hermite interpolation
    Ic(0) = 1: Ic(1) = 1
    Vc(0) = 1 / X(0): Vc(1) = 1 / X(3)
    Call Pchic(Ic(), Vc(), 0, N, X(), Y(), D(), Info)
    If Info <> 0 Then
        Debug.Print "Error in Pchic: Info =", Info
        Exit Sub
    End If
    '-- Compute interpolated values
    Xe(0) = 0.115: Xe(1) = 0.125
    Call Pchfd(N, X(), Y(), D(), Ne, Xe(), Ye(), Yep(), Info)
    Debug.Print "ln(" + CStr(Xe(0)) + ") =", Ye(0)
    Debug.Print "ln(" + CStr(Xe(1)) + ") =", Ye(1)
    Debug.Print "ln'(" + CStr(Xe(0)) + ") =", Yep(0)
    Debug.Print "ln'(" + CStr(Xe(1)) + ") =", Yep(1)
    Debug.Print "Info =", Info
End Sub

Example Results: ln(0.115) = -2.16285581089825

ln(0.125) = -2.07943944206601

ln'(0.115) = 8.6907851599008

ln'(0.125) = 8.00673456704872

Info = 0