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 (xi, fi) 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. 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, 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
