cub_spline_interp_e_cnd
Computes a cubic spline interpolant, specifying various endpoint conditions. The default interpolant satisfies the “not-a-knot” condition.
Synopsis
#include <imsl.h>
Imsl_f_ppoly *imsl_f_cub_spline_interp_e_cnd (int ndata, float xdata[], float fdata[], …, 0)
The type Imsl_d_ppoly function is imsl_d_cub_spline_interp_e_cnd.
Required Arguments
int ndata (Input)
Number of data points.
float xdata[] (Input)
Array with ndata components containing the abscissas of the interpolation problem.
float fdata[] (Input)
Array with ndata components containing the ordinates for the interpolation problem.
Return Value
A pointer to the structure that represents the cubic spline interpolant. If an interpolant cannot be computed, then NULL is returned. To release this space, use imsl_free.
Synopsis with Optional Arguments
#include <imsl.h>
Imsl_f_ppoly *imsl_f_cub_spline_interp_e_cnd (int ndata, float xdata[], float fdata[],
IMSL_LEFT, int ileft, float left,
IMSL_RIGHT, int iright, float right,
IMSL_PERIODIC,
0)
Optional Arguments
IMSL_LEFT, int ileft, float left (Input)
Set the value for the first or second derivative of the interpolant at the left endpoint. If ileft = i, then the interpolant s satisfies
s(i)(xL) = left
where xL is the leftmost abscissa. The only valid values for ileft are 1 or 2.
IMSL_RIGHT, int iright, float right (Input)
Set the value for the first or second derivative of the interpolant at the right endpoint. If iright = i, then the interpolant s satisfies
s(i)(xR)= right
where xR is the rightmost abscissa. The only valid values for iright are 1 or 2.
IMSL_PERIODIC
Compute the C2 periodic interpolant to the data. That is, we require
s(i)(xL) = s(i)(xR) i = 0, 1, 2
where s, xL, and xR are defined above.
Description
The function imsl_f_cub_spline_interp_e_cnd computes a C2 cubic spline interpolant to a set of data points (xi, fi) for i = 0, …, ndata − 1 = n. The breakpoints of the spline are the abscissas. We emphasize here that for all the univariate interpolation functions, the abscissas need not be sorted. Endpoint conditions are to be selected by the user. The user may specify “not-a-knot” or first derivative or second derivative at each endpoint, or C2 periodicity may be requested (see de Boor 1978, Chapter 4). If no defaults are selected, then the “not-a-knot” spline interpolant is computed. If the IMSL_PERIODIC keyword is selected, then all other keywords are ignored; and a C2 periodic interpolant is computed. In this case, if the fdata values at the left and right endpoints are not the same, then a warning message is issued; and we set the right value equal to the left. If IMSL_LEFT or IMSL_RIGHT are selected (in the absence of IMSL_PERIODIC), then the user has the ability to select the values of the first or second derivative at either endpoint. The default case (when the keyword is not used) is the “not-a-knot” condition on that endpoint. Thus, when no optional arguments are chosen, this function produces the “not-a-knot” interpolant.
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, Chapters 4 and 5).
The return value for this function is a pointer to the structure Imsl_f_ppoly. The calling program must receive this in a pointer Imsl_f_ppoly *ppoly. This structure contains all the information to determine the spline (stored as a piecewise polynomial) that is computed by this function. For example, the following code sequence evaluates this spline at x and returns the value in y
y = imsl_f_cub_spline_value (x, ppoly, 0)
The difference between the default (“not-a-knot”) spline and the interpolating cubic spline, which has first derivative set to 1 at the left end and the second derivative set to −90 at the right end, is illustrated in the following figure.
Figure 1, Two Interpolating Splines
Examples
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. Since we are using the default settings, the interpolant is determined by the “not-a-knot” condition (see de Boor 1978).
#include <imsl.h>
#include <stdio.h>
#include <math.h>
#define NDATA 11
/* Define function */
#define F(x) (float)(sin(15.0*x))
int main()
{
int i;
float fdata[NDATA], xdata[NDATA], x, y;
Imsl_f_ppoly *ppoly;
/* Compute xdata and fdata */
for (i = 0; i < NDATA; i++) {
xdata[i] = (float)i /((float)(NDATA-1));
fdata[i] = F(xdata[i]);
}
/* Compute cubic spline interpolant */
ppoly = imsl_f_cub_spline_interp_e_cnd (NDATA, xdata, fdata, 0);
/* Print results */
printf(" x F(x) Interpolant Error\n");
for (i = 0; i < 2*NDATA-1; i++){
x = (float) i /(float)(2*NDATA-2);
y = imsl_f_cub_spline_value(x,ppoly,0);
printf(" %6.3f %10.3f %10.3f %10.4f\n", x, F(x), y,
fabs(F(x)-y));
}
}
Output
x F(x) Interpolant Error
0.000 0.000 0.000 0.0000
0.050 0.682 0.809 0.1270
0.100 0.997 0.997 0.0000
0.150 0.778 0.723 0.0552
0.200 0.141 0.141 0.0000
0.250 -0.572 -0.549 0.0228
0.300 -0.978 -0.978 0.0000
0.350 -0.859 -0.843 0.0162
0.400 -0.279 -0.279 0.0000
0.450 0.450 0.441 0.0093
0.500 0.938 0.938 0.0000
0.550 0.923 0.903 0.0199
0.600 0.412 0.412 0.0000
0.650 -0.320 -0.315 0.0049
0.700 -0.880 -0.880 0.0000
0.750 -0.968 -0.938 0.0295
0.800 -0.537 -0.537 0.0000
0.850 0.183 0.148 0.0347
0.900 0.804 0.804 0.0000
0.950 0.994 1.086 0.0926
1.000 0.650 0.650 0.0000
Example 2
In this example, 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.
#include <imsl.h>
#include <stdio.h>
#include <math.h>
#define NDATA 11
/* Define function */
#define F(x) (float)(sin(15.0*x))
int main()
{
int i, ileft, iright;
float left, right, x, y, fdata[NDATA], xdata[NDATA];
Imsl_f_ppoly *pp;
/* Compute xdata and fdata */
for (i = 0; i < NDATA; i++) {
xdata[i] = (float)(i)/(NDATA-1);
fdata[i] = F(xdata[i]);
}
/* Specify end conditions */
ileft = 1;
left = 0.0;
iright = 2;
right =-225.0*sin(15.0);
/* Compute cubic spline interpolant */
pp = imsl_f_cub_spline_interp_e_cnd(NDATA, xdata, fdata,
IMSL_LEFT, ileft, left,
IMSL_RIGHT, iright, right,
0);
/* Print results for first half */
/* of interval */
printf(" x F(x) Interpolant Error\n\n");
for (i=0; i<NDATA; i++){
x = (float)(i)/(float)(2*NDATA-2);
y = imsl_f_cub_spline_value(x,pp,0);
printf(" %6.3f %10.3f %10.3f %10.4f\n", x, F(x), y,
fabs(F(x)-y));
}
}
Output
x F(x) Interpolant Error
0.000 0.000 0.000 0.0000
0.050 0.682 0.438 0.2441
0.100 0.997 0.997 0.0000
0.150 0.778 0.822 0.0442
0.200 0.141 0.141 0.0000
0.250 -0.572 -0.575 0.0038
0.300 -0.978 -0.978 0.0000
0.350 -0.859 -0.836 0.0233
0.400 -0.279 -0.279 0.0000
0.450 0.450 0.439 0.0111
0.500 0.938 0.938 0.0000
Example 3
This example computes 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, the exact function values are computed with the values of the spline.
#include <imsl.h>
#include <stdio.h>
#include <math.h>
#define NDATA 11
/* Define function */
#define F(x) (float)(sin(15.0*x))
int main()
{
int i, ileft, iright;
float left, right, x, y, fdata[NDATA],
xdata[NDATA];
Imsl_f_ppoly *pp;
/* Compute xdata and fdata */
for (i = 0; i < NDATA; i++) {
xdata[i] = (float)(i)/(NDATA-1);
fdata[i] = F(xdata[i]);
}
/* Specify end conditions */
ileft = 2;
left = 0.0;
iright = 2;
right = 0.0;
/* Compute cubic spline interpolant */
pp = imsl_f_cub_spline_interp_e_cnd(NDATA, xdata, fdata,
IMSL_LEFT, ileft, left,
IMSL_RIGHT, iright, right,
0);
/* Print results for first half */
/* of interval */
printf(" x F(x) Interpolant Error\n\n");
for (i = 0; i < NDATA; i++){
x = (float)(i)/(float)(2*NDATA-2);
y = imsl_f_cub_spline_value(x,pp,0);
printf(" %6.3f %10.3f %10.3f %10.4f\n", x, F(x), y,
fabs(F(x)-y));
}
}
Output
x F(x) Interpolant Error
0.000 0.000 0.000 0.0000
0.050 0.682 0.667 0.0150
0.100 0.997 0.997 0.0000
0.150 0.778 0.761 0.0172
0.200 0.141 0.141 0.0000
0.250 -0.572 -0.559 0.0126
0.300 -0.978 -0.978 0.0000
0.350 -0.859 -0.840 0.0189
0.400 -0.279 -0.279 0.0000
0.450 0.450 0.440 0.0098
0.500 0.938 0.938 0.0000
Example 4
This example computes the cubic spline interpolant to a function, and imposes the periodic end conditions 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.
#include <imsl.h>
#include <stdio.h>
#include <math.h>
#define NDATA 11
/* Define function*/
#define F(x) (float)(sin(x))
int main()
{
int i;
float x, y, twopi, fdata[NDATA], xdata[NDATA];
Imsl_f_ppoly *pp;
/* Compute xdata and fdata */
twopi = 2.0*imsl_f_constant("pi", 0);
for (i = 0; i < NDATA; i++) {
xdata[i] = twopi*(float)(i)/(NDATA-1);
fdata[i] = F(xdata[i]);
}
fdata[NDATA-1] = fdata[0];
/* Compute periodic cubic spline */
/* interpolant */
pp = imsl_f_cub_spline_interp_e_cnd(NDATA, xdata, fdata,
IMSL_PERIODIC,
0);
/* Print results for first half */
/* of interval */
printf(" x F(x) Interpolant Error\n\n");
for (i = 0; i < NDATA; i++){
x = (twopi/20.)*i;
y = imsl_f_cub_spline_value(x, pp, 0);
printf(" %6.3f %10.3f %10.3f %10.4f\n",x,F(x), y,
fabs(F(x)-y));
}
}
Output
x F(x) Interpolant Error
0.000 0.000 0.000 0.0000
0.314 0.309 0.309 0.0001
0.628 0.588 0.588 0.0000
0.942 0.809 0.809 0.0004
1.257 0.951 0.951 0.0000
1.571 1.000 1.000 0.0004
1.885 0.951 0.951 0.0000
2.199 0.809 0.809 0.0004
2.513 0.588 0.588 0.0000
2.827 0.309 0.309 0.0001
3.142 -0.000 -0.000 0.0000
Warning Errors
|
IMSL_NOT_PERIODIC |
The data is not periodic. The rightmost fdata value is set to the leftmost fdata value. |
Fatal Errors
|
IMSL_DUPLICATE_XDATA_VALUES |
The xdata values must be distinct. |
