CSDEC

Computes the cubic spline interpolant with specified derivative endpoint conditions.

Required Arguments

XDATA — Array of length NDATA containing the data point abscissas. (Input) The data point abscissas must be distinct.

FDATA — Array of length NDATA containing the data point ordinates. (Input)

ILEFT — Type of end condition at the left endpoint. (Input)

ILEFT	Condition
0	“Not-a-knot” condition
1	First derivative specified by DLEFT
2	Second derivative specified by DLEFT

DLEFT — Derivative at left endpoint if ILEFT is equal to 1 or 2. (Input)
If ILEFT = 0, then DLEFT is ignored.

IRIGHT — Type of end condition at the right endpoint. (Input)

IRIGHT Condition

0 “Not-a-knot” condition

1 First derivative specified by DRIGHT

2 Second derivative specified by DRIGHT

DRIGHT — Derivative at right endpoint if IRIGHT is equal to 1 or 2. (Input) If IRIGHT = 0 then DRIGHT is ignored.

BREAK — Array of length NDATA containing the breakpoints for the piecewise cubic representation. (Output)

CSCOEF — Matrix of size 4 by NDATA containing the local coefficients of the cubic pieces. (Output)

Optional Arguments

NDATA — Number of data points. (Input)
Default: NDATA = size (XDATA,1).

FORTRAN 90 Interface

Generic: CALL CSDEC (XDATA, FDATA, ILEFT, DLEFT, IRIGHT, DRIGHT, BREAK, CSCOEF [, …])

Specific: The specific interface names are S_CSDEC and D_CSDEC.

FORTRAN 77 Interface

Single: CALL CSDEC (NDATA, XDATA, FDATA, ILEFT, DLEFT, IRIGHT, DRIGHT, BREAK, CSCOEF)

Double: The double precision name is DCSDEC.

Description

The routine CSDEC computes a C2 cubic spline interpolant to a set of data points (xi, fi) for i = 1, …, NDATA = N. The breakpoints of the spline are the abscissas. Endpoint conditions are to be selected by the user. The user may specify not-a-knot, first derivative, or second derivative at each endpoint (see de Boor 1978, Chapter 4).

If the data (including the endpoint conditions) arise from the values of a smooth (say C4) function f, i.e. fi = f(xi), then the error will behave in a predictable fashion. Let ξ be the breakpoint vector for the above spline interpolant. Then, the maximum absolute error satisfies

where

For more details, see de Boor (1978, Chapter 4 and 5).

Comments

1. Workspace may be explicitly provided, if desired, by use of C2DEC/DC2DEC. The reference is:

CALL C2DEC (NDATA, XDATA, FDATA, ILEFT, DLEFT, IRIGHT, DRIGHT, BREAK, CSCOEF, IWK)

The additional argument is:

IWK — Work array of length NDATA.

2. The cubic spline can be evaluated using CSVAL; its derivative can be evaluated using CSDER.

3. Note that column NDATA of CSCOEF is used as workspace.

Examples

Example 1

In Example 1, a cubic spline interpolant to a function f is computed. The value of the derivative at the left endpoint and the value of the second derivative at the right endpoint are specified. The values of this spline are then compared with the exact function values.

USE CSDEC_INT

USE UMACH_INT

USE CSVAL_INT

IMPLICIT NONE

INTEGER ILEFT, IRIGHT, NDATA

PARAMETER (ILEFT=1, IRIGHT=2, NDATA=11)

INTEGER I, NINTV, NOUT

REAL BREAK(NDATA), COS, CSCOEF(4,NDATA), DLEFT,&

DRIGHT, F, FDATA(NDATA), FLOAT, SIN, X, XDATA(NDATA)

INTRINSIC COS, FLOAT, SIN

! Define function

F(X) = SIN(15.0*X)

! Initialize DLEFT and DRIGHT

DLEFT = 15.0*COS(15.0*0.0)

DRIGHT = -15.0*15.0*SIN(15.0*1.0)

! Set up a grid

DO 10 I=1, NDATA

XDATA(I) = FLOAT(I-1)/FLOAT(NDATA-1)

FDATA(I) = F(XDATA(I))

10 CONTINUE

! Compute cubic spline interpolant

CALL CSDEC (XDATA, FDATA, ILEFT, DLEFT, IRIGHT, &

DRIGHT, BREAK, CSCOEF)

! Get output unit number

CALL UMACH (2, NOUT)

! Write heading

WRITE (NOUT,99999)

99999 FORMAT (13X, 'X', 9X, 'Interpolant', 5X, 'Error')

NINTV = NDATA - 1

! Print the interpolant on a finer grid

DO 20 I=1, 2*NDATA - 1

X = FLOAT(I-1)/FLOAT(2*NDATA-2)

WRITE (NOUT,'(2F15.3,F15.6)') X, CSVAL(X,BREAK,CSCOEF),&

F(X) - CSVAL(X,BREAK,&

CSCOEF)

20 CONTINUE

END

Output

X Interpolant Error

0.000 0.000 0.000000

0.050 0.675 0.006332

0.100 0.997 0.000000

0.150 0.759 0.019485

0.200 0.141 0.000000

0.250 -0.558 -0.013227

0.300 -0.978 0.000000

0.350 -0.840 -0.018765

0.400 -0.279 0.000000

0.450 0.440 0.009859

0.500 0.938 0.000000

0.550 0.902 0.020420

0.600 0.412 0.000000

0.650 -0.312 -0.007301

0.700 -0.880 0.000000

0.750 -0.947 -0.020391

0.800 -0.537 0.000000

0.850 0.182 0.000497

0.900 0.804 0.000000

0.950 0.959 0.035074

1.000 0.650 0.000000

Example 2

In Example 2, we compute the natural cubic spline interpolant to a function f by forcing the second derivative of the interpolant to be zero at both endpoints. As in the previous example, we compare the exact function values with the values of the spline.

USE CSDEC_INT

USE UMACH_INT

IMPLICIT NONE

INTEGER ILEFT, IRIGHT, NDATA, NOUT

PARAMETER (ILEFT=2, IRIGHT=2, NDATA=11)

INTEGER I, NINTV

REAL BREAK(NDATA), CSCOEF(4,NDATA), DLEFT, DRIGHT,&

F, FDATA(NDATA), FLOAT, SIN, X, XDATA(NDATA), CSVAL

INTRINSIC FLOAT, SIN