.p>.CMCH3.DOC!CUB_SPLINE_INTERP_E_CND;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 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 3- 1 Two Interpolating Splines
Examples
Example 1
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))
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))
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))
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 functions, 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))
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.
|
Visual Numerics, Inc. PHONE: 713.784.3131 FAX:713.781.9260 |