CSINT
Computes the cubic spline interpolant with the ‘not-a-knot’ condition.
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)
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)
NDATA must be at least 2.
Default: NDATA = size (XDATA,1).
FORTRAN 90 Interface
Generic: CALL CSINT (XDATA, FDATA, BREAK, CSCOEF [])
Specific: The specific interface names are S_CSINT and D_CSINT.
FORTRAN 77 Interface
Single: CALL CSINT (NDATA, XDATA, FDATA, BREAK, CSCOEF)
Double: The double precision name is DCSINT.
Description
The routine CSINT computes a C2 cubic spline interpolant to a set of data points (xifi) for i = 1, NDATA = N. The breakpoints of the spline are the abscissas. Endpoint conditions are automatically determined by the program. These conditions correspond to the “not-a-knot” condition (see de Boor 1978), which requires that the third derivative of the spline be continuous at the second and next-to-last breakpoint. If N is 2 or 3, then the linear or quadratic interpolating polynomial is computed, respectively.
If the data points arise from the values of a smooth (say C4) function f, i.e. ff(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, pages 55 56).
Comments
1. Workspace may be explicitly provided, if desired, by use of C2INT/DC2INT. The reference is:
CALL C2INT (NDATA, XDATA, FDATA, 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.
Example
In this example, a cubic spline interpolant to a function F is computed. The values of this spline are then compared with the exact function values.
 
USE CSINT_INT
USE UMACH_INT
USE CSVAL_INT
 
IMPLICIT NONE
! Specifications
INTEGER NDATA
PARAMETER (NDATA=11)
!
INTEGER I, NINTV, NOUT
REAL BREAK(NDATA), CSCOEF(4,NDATA), F,&
FDATA(NDATA), FLOAT, SIN, X, XDATA(NDATA)
INTRINSIC FLOAT, SIN
! Define function
F(X) = SIN(15.0*X)
! 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 CSINT (XDATA, FDATA, 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 and the error
! 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.809 -0.127025
0.100 0.997 0.000000
0.150 0.723 0.055214
0.200 0.141 0.000000
0.250 -0.549 -0.022789
0.300 -0.978 0.000000
0.350 -0.843 -0.016246
0.400 -0.279 0.000000
0.450 0.441 0.009348
0.500 0.938 0.000000
0.550 0.903 0.019947
0.600 0.412 0.000000
0.650 -0.315 -0.004895
0.700 -0.880 0.000000
0.750 -0.938 -0.029541
0.800 -0.537 0.000000
0.850 0.148 0.034693
0.900 0.804 0.000000
0.950 1.086 -0.092559
1.000 0.650 0.000000