CSPER
Computes the cubic spline interpolant with periodic boundary 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)
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 4.
Default: NDATA = size (XDATA,1).
FORTRAN 90 Interface
Generic: CALL CSPER (XDATA, FDATA, BREAK, CSCOEF [, …])
Specific: The specific interface names are S_CSPER and D_CSPER.
FORTRAN 77 Interface
Single: CALL CSPER (NDATA, XDATA, FDATA, BREAK, CSCOEF)
Double: The double precision name is DCSPER.
Description
The routine CSPER 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. The program enforces periodic endpoint conditions. This means that the spline s satisfies s(a) = s(b), sʹ(a) = sʹ(b), and s"(a) = s"(b), where a is the leftmost abscissa and b is the rightmost abscissa. If the ordinate values corresponding to a and b are not equal, then a warning message is issued. The ordinate value at b is set equal to the ordinate value at a and the interpolant is computed.
If the data points arise from the values of a smooth (say C4) periodic 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 320-322).
Comments
1. Workspace may be explicitly provided, if desired, by use of C2PER/DC2PER. The reference is:
CALL C2PER (NDATA, XDATA, FDATA, BREAK, CSCOEF, WK, IWK)
The additional arguments are as follows:
WK — Work array of length 6 * NDATA.
IWK — Work array of length NDATA.
2. Informational error
|
Type |
Code |
Description |
|
3 |
1 |
The data set is not periodic, i.e., the function values at the smallest and largest XDATA points are not equal. The value at the smallest XDATA point is used. |
3. The cubic spline can be evaluated using CSVAL and its derivative can be evaluated using CSDER.
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 IMSL_LIBRARIES
IMPLICIT NONE
INTEGER NDATA
PARAMETER (NDATA=11)
!
INTEGER I, NINTV, NOUT
REAL BREAK(NDATA), CSCOEF(4,NDATA), F,&
FDATA(NDATA), FLOAT, H, PI, SIN, X, XDATA(NDATA)
INTRINSIC FLOAT, SIN
!
! Define function
F(X) = SIN(15.0*X)
! Set up a grid
PI = CONST('PI')
H = 2.0*PI/15.0/10.0
DO 10 I=1, NDATA
XDATA(I) = H*FLOAT(I-1)
FDATA(I) = F(XDATA(I))
10 CONTINUE
! Round off will cause FDATA(11) to
! be nonzero; this would produce a
! warning error since FDATA(1) is zero.
! Therefore, the value of FDATA(1) is
! used rather than the value of
! FDATA(11).
FDATA(NDATA) = FDATA(1)
!
! Compute cubic spline interpolant
CALL CSPER (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
H = H/2.0
! Print the interpolant on a finer grid
DO 20 I=1, 2*NDATA - 1
X = H*FLOAT(I-1)
WRITE (NOUT,'(2F15.3,F15.6)') X, CSVAL(X,BREAK,CSCOEF),&
F(X) - CSVAL(X,BREAK,&
CSCOEF)
20 CONTINUE
END
X Interpolant Error
0.000 0.000 0.000000
0.021 0.309 0.000138
0.042 0.588 0.000000
0.063 0.809 0.000362
0.084 0.951 0.000000
0.105 1.000 0.000447
0.126 0.951 0.000000
0.147 0.809 0.000362
0.168 0.588 0.000000
0.188 0.309 0.000138
0.209 0.000 0.000000
0.230 -0.309 -0.000138
0.251 -0.588 0.000000
0.272 -0.809 -0.000362
0.293 -0.951 0.000000
0.314 -1.000 -0.000447
0.335 -0.951 0.000000
0.356 -0.809 -0.000362
0.377 -0.588 0.000000
0.398 -0.309 -0.000138
0.419 0.000 0.000000
